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SUMMARY 


Generally accepted solutions for the problems of hypersonic 
flight appear, at the moment, to be centered around the use of 
blunt bodies to minimize the heat-transfer rates. There are, 
however, several other solutions to the problem, and, as part of 
an exploratory study of these solutions, a detailed examination 
has been made of the flow over blunt bodies equipped with a 
spike. These tests, carried out at a Mach Number of about 14 
in the Princeton Helium Hypersonic Tunnel, have investigated 
the effect of varying spike lengths for flat-faced and hemispher- 
ically-nosed axially symmetric bodies. Detailed pressure dis- 
tributions have been obtained as well as heat-transfer rates. 

These exploratory studies have shown that the use of a spike 
protruding from a hemispherical-nosed cylinder at 7 = 14 de- 
creased the pressure level by an order of magnitude and the heat 
transfer to a fraction of that measured on a hemisphere without a 
spike. The general technique appears to hold considerable prom- 
ise for hypersonic flight. 


SYMBOLS 
A = area of copper sector, ft.? 
c = specific heat of copper, B.t.u./Ib.°F. 
Ca = drag coefficient, Cy (AA) sin 6 
Camar = Maximum drag coefficient, 2Cp,,4, (AA) sin 
Cp = specific heat of helium, B.t.u./Ib.°F. 
Cy = pressure coefficient, [(P — P1)/(y/2)P:M?] 
= Maximum pressure coefficient, [(Pr — Pi) + 
(y/2)P1M,?] 
d = body diameter, in. 
h = heat-transfer coefficient, B.t.u./ft.2 sec. °F 
k = density ratio across a normal shock 
K = thermal conductivity, B.t.u. ft./ft.? sec.°F. 
LL = spike length 
m = mass of copper section, Ibs. 
M = Mach Number 
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Preliminary Investigations of Spiked Bodies 
at Hypersonic Speeds' 


SEYMOUR M. BOGDONOFF* anno IRWIN E. VAS** 


Princeton Unwersity 


Nu = Nusselt Number, hd/K 
P = measured pressure, psia 
Pr = Prandtl Number, cpu/K 
q = heat transferred; B.t.u./ft.* sec. 
r = recovery factor, 7, — 7/T) — 7 
Re = Reynolds Number pus/u 
Rs = radius of curvature of shock, in. 
s = distance along surface measured from the stagnation 
point, in. 
t = time, sec. 
r = temperature, °R. 
u = velocity, ft./sec. 
B = velocity gradient du2/ds};=0 
+ = ratio of specific heats for helium 
6 = angle measured from the initia! Mow direction, deg. 
m7 = viscosity, lb. sec /ft. 
Subscripts 
0 = stagnation condition 
1 = free-stream condition 


condition at the edge of the body boundary layer 
b = body 

r = recovery 

: = total conditions behind the normal shock 


INTRODUCTION 


I general aeronautical engineering practice, a great 
deal of effort is expended in making sure that the 
flow does not separate from the bodies which are being 
studied. Contrary to this usual practice, the work 
here is aimed at the deliberate use of separation as a 
solution to the important problem of heat transfer to 
a body traveling at hypersonic speeds. The primary 
work at hypersonic speeds has been concerned with 
the crucial problem of heat transfer, and one way to 
alleviate this has been to use relatively blunt bodies 
to minimize the maximum heat-transfer rate. It was 
recognized quite early by Allen and others' that one 
solution to the heating problem was to use high-drag 
bodies, dissipating a large part of the kinetic energy of 
the body as heat in the air surrounding it. The de- 
creased heat transfer with increase of radius of curva- 
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Fic.1. Princeton Helium Hypersonic Tunnel, settling chamber, 
nozzle, and test section with window removed. 


ture at the stagnation point was also recognized. As 
a consequence, a large amount of work has been in- 
volved with the solution of the flow over a blunt body 
with its high drag and low final speed, surrounded by 
an envelope of very high temperature air which is 
dissociated and ionized. The calculation of this flow 
field, complete with dissociation, ionization, recom- 
bination, surface catalysis, etc., has not been reduced 
to engineering design procedures for general body 
shapes. 

The general idea for the present investigation came 
from supersonic studies of blunt bodies with spikes. 
Two techniques by which the heat transfer can pos- 
sibly be reduced are involved. First, if the separated 
flow remains laminar, the heat transfer to the solid 
wall in the separated region will be only a fraction of 
that associated with attached flow, and second, if the 
flow can be made to separate ahead of a blunt body, 
the resulting external flow will correspond to that for 
a slender body and the pressures over the body will be 
much reduced. This alteration of the basic flow field 
will also result in considerably lower pressures over the 
latter part of the body and a lower base pressure as 
well, thus providing the possibility of decreasing the 
heat transfer over the entire body just because of a 
considerable reduction of pressure level. 

At supersonic speeds, the formation of separation 
ahead of a body caused by a spike has not seen con- 
siderable use, particularly because the heat transfers 
associated with the fore part of the body were higher 
with the spike than without. At a Mach Number 
of about 2 for which these studies were carried out,” 
transition of the free boundary layer from the spike 
occurred and the resulting heat transfers were those 
associated with turbulent rather than laminar flows. 
For Mach Numbers in excess of 10, there are no data 
available on transition Reynolds Numbers, and, with 
the exception of some very exploratory work in refer- 
ence 3, there is no information on which to base an 
analysis of such separation phenomena. If the sepa- 
rated flows remain laminar, and the subsequent expected 
decrease of pressure level is obtained, very significant 
decreases in the heat transfer to the body will result. 
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In addition, for such types of bedies (even at very high 
Mach Numbers) the temperature in the general flow 
field surrounding the bodies is not sufficiently high to 
cause dissociation and ionization. In fact, these 
problems would be restricted only to a very small re- 
gion at the tip of the spike. With the exception of 
this small region, the flow field around the body would be 
purely fluid mechanical and, therefore, should be quite 
susceptible to the general techniques of analysis known 
in supersonic flows. It is, of course, immediately ob- 
vious that, if such a flow could be established, the drag 
of the body would be considerably lower than that of a 
blunt body. Depending upon the extent of the sepa- 
rated region, variable drag characteristics of the body 
could be obtained. The very important assumption 
made here is that it is possible to make such a spike for 
flight applications. This detailed problem is not 
covered in this report, but it is believed that techniques 
are available to permit this to be done. 

The studies reported here were carried out to examine 
the possibilities of obtaining this type of separated 
flow and to get some idea as to whether the possible 
gains, as compared to a body without separated flow, 
might be sufficient to compensate for a somewhat 
more complicated mechanical structure. Two simple 
bodies were examined: a flat-nosed cylinder and a 
cylinder with a hemispherical nose, both with spikes 
of various lengths extending from the nose. Detailed 
examinations of the pressure distributions and _ heat- 
transfer rates have been made on these bodies and 
compared with those for the bodies without spikes. 
The tests were purely exploratory in character and 
were carried out only at zero angle of attack. No at- 
tempt was made to optimize body shape for minimum 
heat transfer or drag, but rather to explore the possi- 
bilities of the general phenomenon. The experimental 
studies were carried out in the Princeton Helium 
Hypersonic Tunnel at Mach Numbers of 12.7 and 
14.0. 


EXPERIMENTAL FACILITIES AND MODELS 


The experimental investigation was conducted in 
the Princeton Helium Hypersonic Tunnel.‘ The 
nozzle consisted of a 5'/2° conical tapered section 
joined to a 3!/,in. diameter cylindrical section. The 
sides of the nozzle were milled parallel to each other 
and windows installed. A photograph of the tunnel is 
shown in Fig. 1. Because of the conical design of the 
nozzle, a Mach Number gradient of about 0.5 per in. 
exists along the center line of the nozzle. The ex- 
periments were conducted at Mach Numbers between 
12 and 14 and at a stagnation pressure of about 1,000 
psia, using helium as the test fluid. 

Two geometrical configurations were studied—a flat- 
nosed cylinder and a cylinder with a hemispherical 
nose. The models were supported in the tunnel by a 
vertical strut attached to the cylindrical portion of the 
model, approximately 6 body diameters from the nose. 
The pressure distribution models were made of brass 
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and were 5/S in. and 1/2 in. in diameter. For tests 
with a spike, the stagnation pressure orifice was elimi- 
nated and a 0.050-in. diameter steel rod with an 11° 
total angle conical nose passed through the body. Its 
position relative to the body was adjusted by a motor- 
driven gear arrangement which permitted the length 
of the spike extending ahead of the body to be varied. 
A photograph of the pressure distribution models is 
shown in Fig. 2(a). The pressures over the models 
were measured by either absolute oil or mercury ma- 
nometers and were recorded photographically. 


Heat-transfer measurements were made by tie 
transient thin wall technique. The temperature of a 
known copper section of the model, as well as the 
stagnation temperature of the tunnel, were measured 
and recorded on a high-speed potentiometer-printc ut 
system. With the available equipment, a record of 
each temperature could be made every 1'/2 sec. For 
the flat-nosed cylinder, a hollow 3/8-in. nylon body 
A copper slug 1/4 in. in diameter and 
A copper 


was employed. 
1/8 in. thick was located in the flat nose. 
constantan thermocouple junction was embedded in 
this plug. Heat-transfer measurements on the hemu- 
sphere cylinder were made with 1 /2-in. diameter models. 
Two models were used for this purpose. For one of 
the models, the overall heat transfer to the entire 
hemispherical section was measured by using a nose 
section of copper 3/32 in. thick. The cylindrical 
section of this model was a hollow nylon tube. The 
second model had a copper sector which covered a 
total angle of 90° of the hemispherical nose. This 
sector is here referred to as the fore-sector. The re- 
mainder of the hemispherical section and the following 
cylinder were of nylon. A photograph of the heat- 
transfer models is shown in Fig. 2(b). 


A conventional two-mirror schlieren system was em- 
ployed for the photographic studies. A spark light 
source of approximately 1 microsec. duration was used 
for most tests, but for special nonsteady studies, a 
continuous light source was used with a Fastax motion 
picture camera operating at 8,000 frames per sec. 


INCHES 


Fic. 2 (a). Photograph of the pressure distribution models. 
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TEST PROGRAM 


The test program was set up to make three types of 
studies: (a) optical, (b) pressure distribution, and 
(c) heat transfer. Detailed pressure distributions were 
obtained over the nose and cylindrical sections of both 
models for a range of spike lengths varying from 0 to 
about 8 body diameters. For each test, spark schlieren 
photographs were taken, and, when unsteady flow was 
suspected, a series of spark photographs and a Fastax 
motion picture were made. For the same range of 
spike lengths used above, heat transfers to the nose 
section were measured using a transient temperature 
technique. 


RESULTS 


(a) Optical Studies 


For each body (diameter, d), spark schlieren photo- 
graphs were taken as the length of the spike, L, was 
varied over the range 0 < L/d < 8. Some selected 
photographs are presented in Figs. 3 and 4 for a free- 
stream Mach Number of 14 at the nose and a Reynolds 
Number of 0.73 X 10° per in. 


(b) Pressure Distribution 


The pressures over the two models were measured 
for the same range of spike lengths investigated in sec- 
tion (a). The results are presented as the pressure 
coefficient ratio C,/Cprmar versus the dimensionless 
distance s/d. 

In Figs. 5(a) and 5(b), the pressure coefficient on the 
nose of the flat-nosed cylinder is shown for various 
L/d. For L/d = 6, the pressure coefficient on the 
nose has reached a point where, for larger values of 
L/d, further decreases in the pressure coefficient are 
negligible. A summary of the nose pressure distribu- 
tion for a few values of L ‘dis shown in Fig. 5(c). The 
pressure on the cylindrical section of the same model 
is shown in Fig. 6 at a Mach Number of 12.7 for values 
of O< L/d <4. ForL/d > 2, the effect of increasing 
the spike length was negligible. 
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Fic. 2 (b). Photograph of the heat-transfer models. 
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Results of the studies on the nose section of the 
hemisphere cylinder are presented in Figs. 7(a) and 
7(b). A summary of the nose pressure distributions 
is shown in Fig. 7(c) for a few selected values of L/d. 
The details of the flow on the cylindrical section of this 
model are shown for 0 < L/d < 5 in Fig. 8. For 
values of 4 < L/d < 9 no change of the pressure 
coefficient occurred. 


(c) Heat-Transfer Measurements 


‘ For the range of spike lengths studied in sections 
: (a) and (b), gross heat-transfer measurements were 
made over the nose at a Mach Number of 14. This 
was done by using copper sections, as described pre- 


Fic. 3. Schlieren photographs of the flow over a flat-nosed 
cylinder with spikes of various lengths; M = 14, Re/in. = 0.73 
106. 
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viously, and measuring their temperature change with 
time in response to a step function imposed by the 
tunnel start. Both preheated and precooled models 
were examined as a check on the heat-transfer results. 
From the temperature measurements, the recovery 
factor and heat-transfer coefficient can be obtained. 
The amount of heat from the copper section is 
g(t) = mc(dT, /dt) 
The amount of heat transferred to the section is 
= hA(T, — T)) 
The heat-transfer coefficient is now defined as 


h (mc/A) [1/(T, — T,)] (dT, dt) 


L/d=7 


Fic. 4. Schlieren photographs of the flow over a hemisphere 


cylinder with spikes of various lengths; M/ = 14, Re/in. = 0.73 
xX 108. 
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Fic. 5(a). 


C,/Cpmaz versus s/d on the nose of the flat-nosed 
cylinder 0 < L/d < 2. 
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Fic. 5(c). 


Summary of the pressure coefficient ratio on the nose of 
the flat-nosed cylinder for various L/d. 


Fie. 6. 


Pressure coefficient ratio on the cylindrical section of the 
flat-nosed cylinder for various L/d. 
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Fic. 7(b). 
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Cy/Comar Versus s/d on the nose of the hemisphere 
cylinder for 4 < L/d < 9. 
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In this equation the heat-transfer coefficient and re- 
covery temperature are unknown. As the heat-transfer 
coefficient is a constant (for all values of ¢), the right- 
hand side of the equation may be solved for various 
values of ¢ such that / is unchanged. The amount of 
energy being transferred from the copper section to 
the stream is obtained from the temperature-time 
gradient at a particular value of ¢. A value of r is 
assumed and 7, calculated. The heat-transfer coelli- 
cient, /, can now be obtained. This process is iterated 
for various assumed values of y and measured dT), /dt 
until a constant value of r(t) is obtained which results 
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the hemisphere cylinder for various L/d. 
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re- in a constant value of h. Because of the desired accu- we 
fer racy of the above calculations, data for large d7), dt and pipes: ig 
ht- T | (the initial part of the time history) were REAR SECTOR 
us used. A summary of the heat-transfer coefficients and 
of recovery factors for both models is given in Tables 1, . \ 
to 2, and 3. A graph of heat-transfer coefficient on the , Ay 
me hemispherical nose as a function of L/d is shown in * . 7 
is Fig. 9. \ 
| \\ 
fed DISCUSSION 
1 Optical Studies 
Its 
A general picture of the flow phenomena can be seen : 
from the schlieren photographs. For Z,d = 0, the 


normal detached shock stood at a distance of about 
0.29d from the flat-nosed body. For 1 < Ld < 3, 
the shock was unsteady, and several spark schlieren 
photographs were taken for each L/d as well as a 
Fastax motion picture at speeds up to 8,000 frames 
| per sec. For L/d = 1, the schlieren photographs 
showed that a strong shock oscillated between the 
spike tip and the body. In similar photographs for 
1!» < Ld < 3, the shock seemed attached to the tip 
of the spike while the downstream sections oscillated, 
| appearing convex or concave. In contrast, the high- 
speed motion picture showed no oscillation of the shock 
pattern, indicating that the phenomenon was occurring 
at frequencies above 10‘. For a spike length about 


| 
4 TABLE 1 
Recovery Factor and Heat-Transfer Coefficient for the Flat- 
Nosed Cylinder at Various L/d Values 
L/d r h 
0 0.982 0.0773 
4 2 0.905 0.0537 
6 0.864 0.0337 
TABLE 2 
Recovery Factor and Heat-Transfer Coefficient for the Total 
Hemisphere at Various L/d Values 
Unheated —_—~PPrreheated Precooled 
L/d r h r h r h 
ad 0 0.929 0.0587 0.927 0.0612 
1 0.846 0.06438 0.853 0.0638 
3 0.851 0.0397 0.850 0.0485 0.852 0.0360 
3 0.859 0.0249 0.863 0.0242 0.915 0.0250 
_ 4 0.863 0.0197 0.864 0.0181 0.850 0.0232 
5 0.863 0.0193 0.870 0.0185 
6 0.867 0.0178 0.871 0.0168 
7 0.865 0.0173 0.875 0.0160 
7 0.870 0.0157 0.875 0.0152 
8 0.870 0.0145 0.880 0.01471 
TABLE 3 
Recovery Factor and Heat-Transfer Coefficient for the Fore 
Sector of the Hemisphere 
Unheated Preheated Precooled 
L/d r h r h r h 
| 0 0.976 0.110 
1 0.943 0.02438 0.895 0.0211 
0.953 0.0244 
2 0.897 0.0176 0.824 0.0214 
3 0.874 0.0127 0.880 0.0111 0.850 0.0109 
J 4 0.870 0.0109 0.870 0.00964 0.860 0.00992 
$ 5 0.870 0.0105 0.864 0.00846 0.845 0.0112 
f 6 0.870 0.00913 0.870 0.00826 0.880 0.00743 


0.870 0.00850 0.860 0.00681 0.850 0.0101 


Fic. 9. Heat-transfer coefficient for various L/d. 


4 times the body diameter, the shock was steady and 
the separated region appeared to start from the cone 
tip. For L/d > 4, the separated flow began on the 
cylindrical section of the spike, forming an angle which 
varied from 8° to 6!/,° as L/d varied from 4 to 8. 
Mair carried out a similar study at a Mach Number 
of 1.96 (see reference 5). With spike extensions of 
0.7 < L/d < 1.3, an oscillation of the flow configuration 
about the body took place at about 6,000 cps. For 
L/d ~ 1.5 the flow separated at the spike tip and the 
oscillation ended. At this condition, a 15° half angle 
conical separated region was formed and the minimum 
drag was measured. Longer spike extensions resulted 
in an increased conical angle of separation and higher 
drag. Similar work done by Beastall and Turner, at 
Mach Numbers between 1.5 and 1.8, indicated that 
the oscillation may be diminished by changing the 
spike tip detail,® but no attempt to study this effect 
was made in the present investigation. 

With the hemispherical nose, for L/d = 0.2 the 
shock about the body seemed to be unaffected by the 
spike. A conical shock was formed at the spike tip 
and moved with the spike tip for 0.56 < L/d < 3. 
When the spike was extended beyond about 4 diameters 
the flow separated from the cylindrical section of the 
spike. Reattachment of the separated flow region 
occurred on the hemispherical section of the body. No 
unsteadiness was noticed for any of the spike lengths. 
The separation angle was measured from the schlieren 
photographs and varied from 6'/,° to 5° for spike 
lengths varying from 4 to 8 body diameters. At low 
Mach Numbers, Mair® and Moeckel’ have shown that 
a conical dead air region, starting at the cone tip, was 
formed for a spike extending about 1.5d from the body. 
For L/d > 1.5, the separated region did not change 
considerably, and the flow separation angle was found 
to be about 17°. 


Pressure Distributions 


For the flat-nosed cylinder, the pressure over the 
nose was nearly constant and equal to the pitot pres- 
sure. The pressure coefficient dropped rather rapidly 
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Fic. 10. Drag coefficient ratio for the spikes of various lengths. 


between L/d = 1 and L/d = 1'/s. It should be 
recalled that, for 1 < L/d < 3, the flow phenomenon 
is unsteady and the pressure distribution cannot be 
considered accurate. For L/d > 6, the change in 
spike length does not change the pressure significantly. 

The pressures on the cylindrical section of a flat-nosed 
cylinder have been presented previously over a range 
of hypersonic Mach Numbers.’ The introduction of a 
spike changed the flow configuration such that there 
was a much smaller expansion around the shoulder. 
This results in a continuously decreasing pressure along 
the cylindrical section (see Fig. 6). Increasing the 
spike length from 2 to 4 body diameters caused only a 
very slight change in the pressure distribution. The 
pressure on the cylindrical section could be calculated 
by considering an equivalent cone formed by the sep- 
arated region, followed by an expansion to the axia! 
direction (along the cylindrical surface). The pressure 
coefficient ratio, as obtained by this method, is slightly 
below ambient in comparison to the measured value of 
about 0.007. 

The measured pressures on the face of the body were 
integrated to form a dimensionless drag coefficient, 
Cy. Asareference, a maximum drag coefficient, Cano; 
was obtained by considering the pressure behind a 
normal shock and the nose area. For the flat-nosed 
body, Ci/Camaz is nearly unity. It is noticed (see 
Fig. 10) that, as L/d varies from 4 to 6, the drag varies 
from 0.03 to 0.02 that of the flat-nosed body. 
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Using the separation angles from the schlieren photo- 
graphs, the pressure on the face in the separated region 
was calculated from the equivalent cone by the method 
suggested by Lees.’ In the steady flow region, there 
is good agreement between the measured and calcu- 
lated pressures. Similar studies at Mach Numbers 
between 1.5 and 2 (see references 5 and 6) gave quite 
different results. The nose drag was reduced to only 
1/3 of its value without a spike for an L/d of the order 
of 1}, 4. 

The pressure distributions over the basic hemi- 
spherical body were similar to those obtained previ- 
ously in reference 8. As the spike was extended, a 
separated region formed and reattachment occurred on 
the hemisphere. The reattachment point moved from 
the 45° station for L/d = 0.5 to 70° for the maximum 
L/d tested. The separated flow caused a strong de- 
crease in pressure over the fore part of the hemisphere. 
Although a peak in the pressure distribution was found 
near the reattachment point, the pressure level over the 
rear of the hemispherical nose and the following cylin- 
drical section was considerably below that found for 
the basic body without a spike. 

The measured pressures were integrated over the en- 
tire hemispherical section and a drag coefficient ob- 
tained. For a Newtonian pressure distribution on the 
basic hemispherical body Ca/Cam,; = 0.5. Similar 
calculations were made with various spike lengths, and 
the variation of nose form drag is shown in Fig. 10. A 
minimum value for the drag coefficient ratio was 
about 0.03 at an L/d of about 6. 

Several previous studies have obtained the pressure 
distribution on a hemisphere cylinder over a wide 
Mach Number range in air and helium. For the fore 
part of the hemisphere (up to about an angle of 75°) 
the pressures are reasonably well predicted by New- 
tonian flow corrected for normal shock losses 


None of the available analytical treatments accu- 
rately predict the pressures on the cylindrical section 
just after the nose. 

The drag on a hemispherical fore body with a spike 
has been measured in the supersonic range both in wind 
tunnels and in free flight (see references 2, 5, 7, and 
10-13). The effect of a spike length of 1'/.d was to 
reduce the drag for Mach Numbers between 1.7 and 
2.7. Longer spikes resulted in increased drag. 


= cos? 0 


Pmax 


Heat Transfer 


Several experimental investigations have been con- 
ducted on hemispherical-nosed bodies for both iso- 
thermal and nonisothermal conditions.'‘~'* Most of 
the theoretical treatments have been restricted to the 
stagnation point region. In this area, reasonable 
agreement is obtained between the experiment and the 
theory of Sibulkin.'’ For incompressible laminar flow, 
the heat transfer at the stagnation point is 


Nua = 0.763(Bd?/v)°*> 
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SPIKED BODIES AT HYPERSONIC SPEEDS 


where @ is the velocity gradient di2/ds},.o. The deter- 
mination of 8 is of critical importance in the calcula- 
tion of 

The potential flow solution on the hemisphere for 8 
is 3u;/d. Two solutions for the real flow are those of 
Li and Geiger, '* 

B = (u/R,) V2[1 — (1/2k)]/k 
and Probstein,!® 


B = [w/(d/2)] V2[1 — (1/28) 


The only difference between these two solutions is the 
use of the nose radius by Probstein in place of the shock 
radius used by Li and Geiger. 

In previous studies,* R, was measured in the present 
range of Mach Numbers and found to be equal to 1.4d/2. 
Using this value, the dimensionless velocity gradient 
8 d/u, at the stagnation point is 0.953 (Li and Geiger) 
and 1.332 (Probstein). From the pressure distribu- 
tion measurements of the present work, the velocity 
gradient at the stagnation point can be obtained (Fig. 
11). The resulting magnitude of 6 d/u is equal to 
1.29, which is close to Probstein’s prediction. 

Using the value of # measured on the fore sector, 
conditions behind the shock, and the measured velocity 
gradient, the value of Nu_/(8d?/v)°-’ was equal to 0.51. 
This is about 20 per cent lower than the Sibulkin stag- 
nation-point value, and, considering the extent of the 
heat transfer surface, indicates a reasonable check of the 
theory and experiment. 

For the flat-nosed cylinder, the potential flow solu- 
tion gives a velocity gradient, 6, equal to 4u)/drz. 
This gradient has also been calculated by Probstein,'® 
who obtained 


= (4/34) V2[1 — (1/2k)]/k 


The pressure measurements on the flat-nosed cylinder 
were not sufficiently detailed to permit a good evalu- 
ation of 8. However, using the Probstein value of 6 
with the measured value of 4 gave a value of the di- 
mensionless heat-transfer parameter 
equal to about 0.47. 

Previous studies of heat transfer to blunt bodies with 
a spike have only been carried out at relatively low 
Mach Numbers. For J/ equal to 2.67 and a Reynolds 
Number of 0.285 X 10° (based on free-stream condi- 
tions and body diameter’), it was found that the 
heat-transfer coefficient on the spiked body was about 
twice that for the basic body. The present work shows 
quite a different trend at 17 = 14 and R, = 0.36 X 
10°. 

The heat-transfer coefficient was measured over the 
fore sector and over the entire hemispherical nose for 
various spike lengths with preheated, precooled, and 
unheated initial conditions. For the hemisphere with 
a spike of L/d = 4, the value of 4 drops to about one 
third of the value found without a spike. On the 
fore sector of the hemisphere, at the same value of 
L/d, h decreases to about one tenth of the unspiked 
value. Since # was obtained for the total hemisphere 
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Fic. 11. Velocity over the hemisphere at M = 14. 


and fore sector separately, it was possible to estimate 
an average value of / over the rear sector of the hemi- 
sphere, as shown in Fig. 9. For optimum spike 
lengths, the local heat transfer is lower than for the 
body without a spike. 

Somewhat similar results were found for the nose of 
the flat-nosed cylinder, although only limited tests 
were made. For L/d = 6, h was reduced to less than 
one half the stagnation value. 


General 


The preliminary work that has been detailed in the 
present report indicates that certain advantages can 
be obtained by using a spike protruding ahead of the 
body. For the flat-nosed cylinder, the fore pressure 
drag can be decreased to about one fortieth and the 
heat-transfer coefficient to about one half its stag- 
nation value. The measurements on the hemisphere 
cylinder indicate that the fore pressure drag and aver- 
age heat-transfer coefficient over the entire hemisphere 
can be decreased to about one tenth and one third, 
respectively, of the unspiked body case. The sharp 
decrease in heat transfer with a spike, in contrast to 
the experience at low Mach Number, indicates the 
separated flow is undoubtedly laminar (although no 
direct evidence of this was obtained in the present 
study). These results and results of reference 3 on a 
cone appear to give considerably higher limits of lami- 
nar transition than are expected at lower Mach Num- 
bers, although results at Mach Numbers from 5 to 9 
show this definite trend. 

Several comments are perhaps in order concerning 
the use of the information presented here. The com- 
parison of a spiked body with an unspiked body de- 
pends, in some respects, on the type of application. 
The work reported here shows that a heating rate can 
be reduced concurrently with a reduction in the drag 
coefficient. For a re-entry body, the reduction in heat- 
transfer rate is just about compensated for by the de- 
crease in drag coefficient resulting in approximately the 
same heat input for the spiked and unspiked body. 
In this case, the main advantage to be obtained by the 
use of a spike is that, for a given heat capacity of the 
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body, re-entry can be accomplished at a considerably 
higher Mach Number. This may be of some impor- 
tance in the problem of dispersion and counter measures. 
On the other hand, if one considers a body flying in the 
atmosphere—e.g., a boost glide vehicle which travels 
at approximately a constant Mach Number—the use 
of a spike has a quite different effect. It results in a 
quite considerable decrease in the drag of the body, 
thus perhaps providing an increased L/d ratio or less 
thrust needed to sustain flight, and, for flight at a 
constant Mach Number, results in a reduced heat 
transfer. There is also the possibility that, in some 
cases at least, the application of a variable drag coeffi- 
cient with a body of essentially constant weight and 
frontal area would have several advantages. Such dis- 
cussions have, in the past, been based on the use of 
skirts. It is suggested, however, that the use of a 
variable spike may be considerably simpler than a 
variable geometry skirt to obtain a variable drag co- 
efficient in flight. These examples are used simply 
to point out that the use of a spiked body may be quite 
wide with quite different end results obtained from its 
use for different types of applications. 

Although the heat-transfer measurements are limited 
to the nose section, it is expected that further studies 
will also show lower heat transfers over the cylindrical 
section and base of the spiked body. For spike lengths 
of about L/d = 4 and higher, the entire pressure level 
over the cylindrical section was considerably lower 
than for the unspiked body. The high-pressure levels 
over bodies with blunt noses have been studied in de- 
tail in references 8 and 20, where it has been shown 
that the high-pressure level will be experienced many 
nose diameters downstream on the body following. 
‘The use of the spike simply localizes the region of high 
pressure to the tip of the spike, permitting an order of 
magnitude decrease in the general pressure level around 
the body. With a lower pressure ahead of the base, 
lower base pressures would be expected?! with a prob- 
able lowering of the heat transfer in this region. 


CONCLUSIONS 


On the basis of the series of experimental studies 
which have been carried out on blunt bodies with a 
spike, several important conclusions have been reached. 
It should be noted that the results were obtained in 
helium (‘‘perfect’’ gas with viscosity) and the direct 
problems of temperature associated with the test 
Mach Number of 14 are not included. 

(1) The use of a spike protruding from a hemi- 
spherical-nosed cylinder at .1/ = 14 decreased the 


pressure level by an order of magnitude and the heat 
transfer to a fraction of that measured on a hemisphere 
without a spike. 

(2) For the hemisphere cylinder with optimum spike 
length, L/d ~4, the form drag was one tenth and the 
heat transfer one third of the basic values without a 
spike. 
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(3) With the flat-nosed cylinder and L/d ~4, the 
form drag was one fortieth and the heat transfer one 
half of the section without a spike. 
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A Linearized Analysis of the Forces 
Exerted on a Rigid Wing by a Shock Wave 


F. EDWARD EHLERS anp E. M. SHOEMAKER 
Boeing Airplane Company* 


SUMMARY 


Solutions are obtained in closed form for the pressures exerted 
on a rigid half plane by an incident, plane acoustic shock wave. 
The angle of incidence of the wave front is arbitrary and the half 
plane is considered to be traveling at constant velocity, subsonic 
or supersonic with respect to the acoustic medium. A closed- 
form solution is obtained also for a rigid wedge which is motion- 
less with respect to the acoustic medium. The analysis is carried 
out by transforming the wave equation to Laplace’s equation 
by the Busemann conical transformation and then applying con- 
formal mapping. 


SYMBOLS 


c = velocity of sound 
“a = pressure coefficient 
M = Mach Number 
as = pressure 
U = velocity of wave front 
Uo = perturbation velocity potential in regions (3) 
and (4) of Fig. 6 
u = perturbation velocity behind shock front in 
region (1) 
ee = Cartesian coordinates in physical plane 
x = X/et 
= conical coordinates 
r,@ = cylindrical coordinates in x, y plane 
£ = coordinate axis in Figs. 8 and 9 
(E, n) = final plane and coordinates after Busemann 
(p,w) transformation and conformal mappings 
po = air density 
we 
wp = angles defined in Fig. 2 
wp* 
wQ 
Subscripts 
i refers to the Busemann plane 
2 refers to plane after conformal transformation 


from Busemann plane 


INTRODUCTION 


i hen CONSIDERATION of the loads exerted on an air- 
craft by the blast wave from an exploded bomb is 
important from the points of view of both structural 
design and general weapons analysis. It is extremely 
difficult to obtain accurate experimental data, and 
efforts in this direction have not yielded satisfactory 
results. Therefore, a mathematical theory is highly 
desirable. 

The following analysis is based on linearized theory, 
and all solutions are obtained in closed form. (The 
exact physical problem, involving nonlinear equations 
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and undefined domains of influence, is sufficiently com- 
plicated that an analytic solution is probably not feasi- 
ble.) The linearized theory assumes a weak shock 
such that the velocity and pressure fields behind the 
shock front are small perturbations from the state of 
rest. Both the velocity potential and pressure then 
satisfy the acoustical wave equation. 

In the following, we find the solution of the pressures 
produced by a weak shock front striking a rigid half 
plane which may be regarded as an approximation to 
an aircraft wing. The incidence angle of the shock is 
arbitrary, but the shock front is assumed to be always 
parallel to the edge of the half plane so that the problem 
is two dimensional. The half plane is considered to 
be moving at constant speed in its own plane, sub- 
sonically or supersonically with respect to the undis- 
turbed air. An analysis is made also of a plane shock 
striking a rigid wedge which is at rest with respect 
to the undisturbed air. 

The consideration of a half plane instead of a plate of 
finite chord is not particularly restrictive since, in the 
actual case of a wing with finite chord, the conditions 
of most practical importance in design—i.e., maximum 
chordwise and spanwise bending—should not occur 
after the shock front passes the wing. 

Graphical results are presented of the distribution of 
net pressure (lift) for various shock incidence angles 
and various wing Mach Numbers. 


THE SHOCK WAVE STRIKING A HALF PLANE TRAVELING 
SUBSONICALLY 


Fig. 1 illustrates the geometry of the wave fronts 
in the x, y plane where the variables Y, Y in the physical 
plane have been transformed by the conical coordinate 


transformation x = X/ct, y = Y/ct. Since all dis- 
F 
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Fic. 1. The wave pattern in conical coordinates for a plane 
shock wave striking a half plane traveling subsonically. 


|_| 
_ 
|| 
if a 
mic 
the 4 
luff 
nic 
Re- 
unt 
the 
mi- 
ta 
of 
\S, 
: 
che 
D2. 
dy 
ol. 
Re- 
ds, 
ng 
ds, 
: 
H 


Fic. 2. The definition of the angles wp, wp*, and we for a general 
point w. 


turbances travel with constant velocity, the geometry 
of Fig. 1 is invariant with time and is unique. The 
radius of the circle is unity and the distance OA’is M, 
where / is the Mach Number of the half plane. In 
Fig. 1, A lying to the right of point O corresponds to 
the half plane moving so that the shock overtakes the 
trailing edge, while A to the left of O corresponds to the 
shock striking the leading edge. 

Regions (0). (1), and (2) are regions of constant state. 
From Huyghen’s principle the total energy in a wave 
front remains constant, although the wave front area 
may vary with time. It follows that across the incident 
wave fronts EF and CG and across the reflected wave 
front CH there are pressure discontinuities [P] = ¢; 
across the cylindrical wave front FBHF the pressure is 
continuous because this wave front originated with 
zero area. For convenience, all pressures are given as 
their differences from the pressure in region (0) and are 
normalized to the magnitude of the pressure jump 
across the shock. Consequently, pressures in regions 
(Q), (1), and (2) are equal to 0, 1, and 2, respectively. 

The pressure field is unknown only in the interior of 
the circle of Fig. 1. The boundary conditions on r = 1; 
namely, 

P=0 for CSO0<6 

P=1 for <0<27-0 

P=2 for 2r—>0, < @ < 2s 
follow from the requirement of pressure continuity 
across the unit circle. On the lines AB, and AB_, the 
condition that the flow velocity be tangent to the half 
plane is satisfied by requiring the normal derivative 
of the pressure to vanish. 


The pressure satisfies the wave equation, which, after 
transformation to conical coordinates, becomes 


(1 — r?)Prr + (1/r — 2r)Pr + (1/r?)Po = 0 (2.1) 


where r = Vx? + y? = V (X/ct)? + (Y/ct)? and 6 = 
tan~'(y/x). (2.1) is transformed to Laplace’s 
equation by applying Busemann’s transformation; 
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r= (x? + + (1 — x? — = 
(1 — (2.2) 
6 = (y/x) 


Note that the angle @ and the boundary of the unit 
circle are invariant. 

The general solution for the pressure field inside the 
unit circle in the Busemann plane, Z; = ne”, now may 
be expressed as the real or imaginary part of an analytic 
function of the complex variable Z,; thus, we are able 
to make effective use of conformal mapping. First, 
the unit circle in the Z; plane is mapped into the unit 
circle in the Z: plane with the real axis mapped into the 
real axis and with the point, 


= M/[1+ (1 — M)'?]=B (2.3) 


corresponding to x = M, transformed into the origin. 
Thus, 


= (4, — B)/(l — (2.4) 


The angles +6; of Fig. 1 then transform into +62, given 
by 


+6. = +0, + [@ sin @,/(1 — Bcos@,)] (2.5) 
The next transformation 
w=Z,'", w=pe? in (2.6) 


maps the unit circle in the Z, plane into the unit semi- 
circle in the upper half plane. The lines AB, and 
AB_ in Fig. 1 map into A’B’, and A’B’_, respectively, 
of Fig. 2. 

For convenience, we consider the problem for the 
entire unit circle in Fig. 2 with the boundary conditions 
reflected about the real axis. (Note that for this par- 
ticular problem, in Fig. 2 w = and = — 42/2. 
The solution will be obtained first for general angles 
w, and ws.) By the principle of superposition P is 
taken to be the sum of a pressure ?; and a pressure P». 
The appropriate boundary conditions are as follows. 
For 


P, = im {f(w)}, | 

p=1, -m 

on Sw < 2r — w 

1, < we 

on | 1, 27 —w <w < — w 

For P2: 
P, = Im {g(w)}, <1 
P,=0 on p=1, -mlw< w 


on p=1, w<w 27 — w 


Because of the symmetry about the real axis of the 
boundary conditions on the unit circle, the condition 
that the normal derivative of the pressure vanish on the 
real axis is assured. 

The solution to P; has been found by Ludloff! and is 
given by 


and f 
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= (1/m) [we — (w2 — w)/2] + 
(1/m) [wp* — (we — w)/2] (2.7) 


where wp and wp* are defined in Fig. 2. It is easily seen 
that this function satisfies the boundary conditions for 
P, since an angle inside a circle with its vertex on the 
boundary is equal in radians to one-half its subtended 
arc. Eq. (2.7) is also a harmonic function since it is 
the imaginary part of the analytic function 


w— ei” w— 


2.8) 


fiw) = (1/m) log | 


On the real axis, we = wp* = wp(é, 0) as shown in 
Fig. 3. Then 


0) = — [sin ‘(cos w — — 
tan—! }sin — w)/[cos — + —]} (2.9) 


or 

| sin(w. — wi) + &(sin — sin we) | 

an—! | — 

—cos (a. — + E(cos w + COS we) — 
(2.10) 

For the particular case considered, w, = (7 — w) = 
so that 


wp(t. 0) = — [sin 62/(cos 6. — (2.11) 
and from Eq. (2.7) the pressure P; on the real axis is 


0) = ] + = 

(2/7) tan~'[ sin 62/(cos 6. — (2.12) 
where @2 is given by Eq. (2.5). Finally, € is expressed 
in terms of the variable x by 


1/2 


= = (x, — B)'”?/(1 — (2.13) 
where 

= + (1 — 

B= (1 — 


which follow from Eqs. (2.6), (2.4), and (2.2). 
It is easily shown that the function 


P» = (2/m) [wg — (x — (2.14) 


where wg is defined in Fig. 2, satisfies the appropriate 
boundary conditions for P:. In addition, P: is the im- 
aginary part of the analytic function 


g(w) = (2/m) X 
log [(w — e')/(w — e')]} (2.15) 
Hence, Eq. (2.14) is the solution to P». 
On the real axis, wo = we(t, 0), as shown in Fig. 3. 
Then 
tan [we(é, 0) /2] = sin — we)/[E + cos (r — w) | 
(2.16) 


For the particular case considered here, where a, = 
(3 we) = 2, 


Fic. 3. The geometry in the w plane for a point lying on the 
axis. 


tan [we(é, 0)/2] = sin (02/2) /[E + cos (2.17) 


and from Eq. (2.14) the pressure P, on the real axis be- 
comes 


0) = (2/m) X 
(2 tan—! {sin (6./2)/[é + cos (@./2)|} — 62/2) (2.18) 
where & is expressed in terms of x by Eq. (2.13) and 6 
is given by Eq. (2.5). 

The net pressure on the half plane is affected only 
by P2 since P; is symmetric about the n axis. Conse- 
quently, 


0) — P2(é, 0) 
(4/m)tan~! [(sin 62/2) /(cos 6/2 — &)] — 
(4/mr)tan~! [(sin 62/2) / (cos + 


AP 


Combining the inverse tangents yields 
AP = (4/m)tan—! [(2¢ sin 6/2) /(1 — (2.19) 


As an example, we compute the net pressure on the 
stationary half plane. Then 6 = 0, @ = 6;, and the 
relation between £ and the variable x is given by 


+ (1 — (2.20) 


The net pressure AP is plotted against x in Fig. 4 for 
values of 6; = 10°, 45°, and 90°. The length L of the 
half plane segment BC of Fig. 1 with net pressure equal 
to 2is 


L = secé@ — 1 (2.21) 


This is seen to vanish for the shock normal to the half 
plane (6; = 0) and to become infinite for the shock paral- 
lel to the half plane (6, = 7/2). 

Differentiating Eq. (2.19) with respect to x and using 
Eq. (2.20) reveals that the slope of the curve for net 
pressure becomes infinite at x = Oandx = 1. This is 
also apparent from Fig. 4. 

The net pressure produced by a shock wave inclined 
at the angle 6; = 7/4 striking the trailing edge of a half 
plane receding at Mach Numbers of 0.2, 0.4, 0.6, and 
0.8 is shown in Fig. 5. The character of the distribu- 
tion differs very little from that of the stationary half 
plane. 

The above formulas were derived on the assumption 
that the half plane is moving with subsonic velocity 
and the shock overtakes the trailing edge. The re- 
sults are readily extended to the problem of the leading 
edge of a half plane traveling subsonically toward the 
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The distribution of net pressure on a stationary half 
plane for various shock incidence angles. 
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Fic. 5. The distribution of net pressure ona half plane traveling 
at various subsonic Mach Numbers for a shock incidence angle of 
45°, 
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Fic.6. The wave pattern in conical coordinates for a plane shock 
wave striking a half plane traveling supersonically. 
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Fic. 7. The distribution of net pressure on a half plane travel- 
ing at various supersonic Mach Numbers for a shock incidence 
angle of 45°. 
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Fic. 8. The wave pattern in conical coordinates for a plane 
shock wave striking a stationary wedge where the wave front is 
reflected only from one face of the wedge. 
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Fic. 9. The wave pattern in conical coordinates for a plane 
shock wave striking a stationary wedge where the wave front is 
reflected from both faces of the wedge. 
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FORCES EXERTED ON A RIGI 


shock wave. The pressures are again given by Eqs. 
(2.12), (2.18), and (2.19) but with 8 in Eq. (2.13) re- 
placed by — 3. 


THE SHOCK WAVE STRIKING A HALF PLANE TRAVELING 
SUPERSONICALLY 


Consider the half plane moving into the shock wave 
with supersonic Mach Number .1/; then the regions of 
influence of the disturbances have the form indicated 
in Fig. 6, where regions (0), (1), (2), (3), and (4) are of 
constant state with the pressures in (0), (1), and (2) 
having the same values as for subsonic Mach Numbers. 
The angle u 1s defined by 


= cos! (1/M) (3.1) 


where .J/ is the Mach Number of the half plane. The 
solution inside the unit circle can be obtained in the 
same manner as the subsonic case if the pressures on the 
arcs B’F’ and B’H’ are known. Because of the pres- 
sure jump at the shock front, there is a perturbation 
velocity in region (1) in the direction normal to the 
wave front EF. The magnitude uw of the perturbation 
velocity is determined from 


€ = pocu (3.2) 


where ¢€ is the jump in pressure across the wave front, 
p, the air density, and c the velocity of sound. The 
half plane may then be considered as traveling in the 
region behind the shock with supersonic velocity at an 
angle of attack. Let U be the velocity of the half plane; 
then this angle of attack is 


a = tan~! sin 0,)/(U + u cos ~ (usin 


In terms of the pressure jump in Eq. (3.2) this may be 
written 


a = (€/ (sin 6) (3.3) 


According to linearized supersonic theory, the flow 


in regions (3) and (4) are uniform. If the perturbation 


0) = — 
cos 


from which 


(x1, 0) = 7 


or, since = (1/J/), 


rPi(m, 0) = + 6, + cos! (1/M) — 2 tan 


The boundary conditions that P» satisfies (referring 
again to the analysis for the subsonic case) are P2 = 0 
on are H’BF’ and P, = [1 — (MV sin 4 VM? — 1)] on 
Then, from Eq. (2.18) with 6/2 replaced 


arc H’B’F’. 
by wu, we have 


w.+ 6 — 2 tan ‘| 
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velocity potential in region (3) is U’¢, then 


is a solution of ¢,, + (1 — ./*)¢,, = 0, which satisfies 
the boundary condition of tangential velocity of the 
wing, 


OV = —a 


and of ¢ = 0on AF’. The resulting pressure coefficient 
is 


C, = —26, = —2a/V M? — 1 (3.5) 
and the normalized pressure in region (3) is 
Po = 1+ pol2C,/2e = 1 — M2 — 1 = 
Msing Vir 1 (3.6) 
Similarly, on the lower side, 
Pw =1+ Msin6,/VM — 1 (3.7) 


With the normalized pressure known on the ares of 
the unit circle, we may solve the problem by first divid- 
ing it into two parts, since the disturbance in the upper 
semicircle is not influenced by the disturbance in the 
lower semicircle. We shall consider first the upper 
semicircle. By comparing Fig. 6 with Fig. 2, it be- 
comes apparent that we can make immediate use of 
the solutions for P; and /,, defined by Eqs. (2.7) and 
(2.14), with the angles w and a» of Fig. 2 given by 


w=r—u (3.8) 


The boundary conditions which P; will satisfy are P, = 
1 on are FF’ and P; = 0 on arcs B,F and B_’F’ of Fig. 
6. 
From Eq. (2.9), 
wp(x1, 0) = 2a — [sin 6,/(cos 0; — x1)| — 
[sin w/(cos uw + (3.9) 
(Note that here the z; or Busemann plane is the same 


as the w plane in the discussion of the subsonic prob- 
lem.) Combining the two inverse tangents yields 


sin (6; + + x (sin — sin | (3.10) 
(A; + w) + x1 (cos 0; — cos uw) — 


sin (0; + w) + x (sin #; — sin p) || (3.11) 


cos (0; + + x, (cos — cos pw) — 


(1 + Mx,) sin 0, + (cos 0; — x1) WM? = 1 


(3.12) 
(1 + (cos 6; — — Vi — 1 sin A; 
rP2(x;, 0) = (1 — M sin M2 — 1) X 
tan-! [sin u + — ul (3.13) 


el- 
| 
e 
iS 
kay, 
ge 


or 


aP(x;, 0) = (1 — M sin 6,/W M2? — 1) x 
{2 — 1/(1 + Mx,)] — 


cos~! (1/M)} (3.14) 


The total normalized pressure on the upper side B,’B, 
is then the sum of Eqs. (3.12) and (3.14). 

The lower semicircle solution may be found by using 
the solutions given for the upper semicircle. On 


M sin 6 
= —26, — 2- di cos”! (1/M) 


1 
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B_H, P = 2;0n HH’, P = 1; and on H’B_’, P = 1 + 
M sin 6, VM?—1. We recognize that 


Power =2- 


(3.15) 


Hence the pressure on the lower side B_’B_ follows 
immediately from Eqs. (3.12) and (3.14). 


The net pressure on the plane is found from 


Combining Eqs. (3.12) and (3.14), we get 


M sin 


1+ Mx, 


(1 + sin + (cos 6 — 1) 


The coordinate x; is related to the conical coordinate 
x by Eq. (2.13). The net pressure distribution is shown 
in Fig. 7 for the half plane traveling at Mach Numbers 
of 1.2, +/2, and 2 into a shock wave inclined at 6; = 
45°. Note that there is little variation with Mach 
Number of the pressure distribution on the half plane 
inside the unit circle. However, the length of the 
region of uniform pressure just behind the leading edge 
and the magnitude of this pressure show appreciable 
variation. 


THE SHOCK WAVE STRIKING A STATIONARY WEDGE 


For the shock wave striking a stationary wedge the 
regions of influence have the form indicated in Fig. 8 
for the case 6; > y, or when the point of tangency of the 
shock front to the unit circle (point F) is exterior to the 
wedge. This solution was found previously by Keller 
and Blank? but is included for the sake of completeness. 
Since in the x, y plane the wedge lies on radial lines of 
the unit circle, the boundary of the unit circle and the 
wedge outline will be unchanged by the Busemann 
transformation. Then the transformation 


w = (4.1) 
takes the portion of the unit circle exterior to the wedge 
into the unit semicircle in the upper half plane as in 
Fig. 2. The boundary conditions are identical with 
the case of the half plane traveling subsonically. It 
need only be noted that in this case the angles a; and 
are 


= 1(6; (27 = (4 2) 


w = — — 


The pressure on the wedge is the sum of pressures 
P, and P:. Since on the half plane wp* = wp, Egs. 
(2.7) and (4.2) yield 


P; = —1 + (20, — y)/(2e — y) + 
(2/m)wp(—, 0) (4.3) 


where wp(é, 0) is given by Eq. (2.10) with w; and w, 
from Eq. (4.2). 


Eqs. (2.14) and (4.2) lead to 


(1 + (cos 6 — x:) — VM? — 1 sin 6 


Py = (2/m) [welt 0) — — y)] (4-4) 


where wo(é, 0) is given by Eq. (2.16) with w, from Eq. 
(4.2). 

The variable é is related to the variables x and < of 
Fig. 8 by 


_ sh(x), «20 
z>0 (4.5) 


where 
h(x) = {x/[1 + (1 — (4.6) 


When 6; < y the shock wave is reflected also from 
the upper surface of the wedge and the regions of in- 
fluence are shown in Fig. 9. Again, the transformation 
(4.1) is used with the same results. The angles w and 
are 


= — — yt 


The boundary conditions that the pressure satisfies 
on the unit semicircle of Fig. 2 are P = 2 on B,’F’, P 
= 1 on F’H’ and P = 2 on H’B_’. We recognize 
that the pressure in the unit circle of Fig. 2 is given by 
P = 2 — P, where P, is given by Eq. (2.7). From 
Eqs. (2.7) and (4.7) the solution is 


P =3— — y) — (2/m)wr(é, 0) (4.8) 


where wp(é, 0) is given by Eq. (2.10) with w. and w; 
from Eq. (4.7). The variable & is related to vari- 
ables x and £ of Fig. 9 by Eqs. (4.5) and (4.6). 


CONCLUSIONS 


Formulas have been derived for the pressure exerted 
on a half plane by an incident weak plane shock. The 
shock front was considered to be parallel to the edge 
of the half plane, but the angle of incidence was arbi- 
trary. Three cases were considered: the half plane 
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On Panel Flutter in the Presence of a 


Boundary Layer' 


JOHN W. MILES* 
University of California at Los Angeles 


SUMMARY 


The energy transfer from the uniform flow outside a boundary 
layer to a transverse surface wave at the boundary is calculated 
on the hypothesis that the boundary layer may be represented by 
an inviscid, approximately parallel shear flow. This energy trans- 
fer is found to consist of two components. The first is similar 
to that found previously for supersonic panel flutter in the ab- 
sence of a boundary layer and is relatively diminished by the 
presence of the boundary layer. The second is intrinsic to the 
shear flow and is present whenever the surface wave speed is 
smaller than the free-stream speed (whether subsonic or super- 
sonic). Approximate solutions to the differential equations for 
the disturbed motion of the boundary layer are established and 
compared with more exact numerical solutions. The application 
of the results to actual flutter calculations is discussed, and an 
example of a pressurized cylinder is considered. It is concluded 
that the presence of a typical boundary layer may reduce the 
degree of instability for supersonic flutter of a long monocoque 
shell by an order of magnitude, thereby providing a possible ex- 
planation of discrepancies between earlier theories (in which 
boundary-layer effects were neglected) and observation. 


(1) INTRODUCTION 


_ DYNAMIC INSTABILITY (flutter) of an infinite 
panel exposed to a supersonic air stream has been 
discussed in references 1 and 2. It was found there 
that the flutter speed was approximately a; + (Vo)min 
(where a; denotes sonic speed and Vo the speed of a 
transverse wave on the free panel), and that the critical 
wavelengths for monocoque structures were sufficiently 
small to justify the neglect of panel curvature on the 
aerodynamic (but not the structural) forces and to 
justify the assumption of an infinite panel. The small- 
ness of the wavelength implies that the effects of a rela- 
tively thick boundary layer may be important, and the 
following analysis is intended to establish the qualita- 
tive nature of boundary layer effects and to provide 
first approximations to their magnitude. 

We adopt as a model for our calculations an inviscid, 
parallel shear flow over an infinite, approximately plane 
panel. Assuming a periodic surface wave displacement 
of the panel, we shall establish the equations of motion 
for this shear flow and the resulting pressure on the 
boundary in Section (2), after which the eigenvalue 
problem for the panel will be established in Section (3). 
We then shall consider, in Section (4), the general 
nature of the possible instability and show that, in addi- 
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tion to the essentially supersonic instability of reference 
2, there also appears a new type of instability associated 
with energy transfer from the shear flow; this latter 
mechanism is closely related to that which enters the 
boundary-layer stability problem. 

The boundary value problem considered in Section 
(2) leads to a second-order differential equation that 
has a regular singularity in the range of integration if 
the surface wave speed is less than the air speed at the 
outer (free-stream) edge of the boundary layer. We 
consider, in Section (5), the nature of the solutions to 
this differential equation in the neighborhood of its 
singularity, its solution by expansion in powers of the 
ratio of boundary-layer thickness to wavelength, and 
its solution by numerical integration. The application 
of these results to actual flutter analyses is discussed in 
Sections (6) and (7). 

We emphasize that the model on which the present 
analysis is based is highly idealized, and experimental 
confirmation of the results is most desirable; moreover, 
even though the predicted conditions for the initiation 
of flutter should prove correct, nonlinear effects must 
be decisive in determining the severity of such flutter. 
In the absence of experimental confirmation and at 
least a qualitative understanding of these nonlinear 
effects (especially with reference to structural fatigue), 
our results must be accepted with considerable reserve 
and should not be regarded as adequate for major de- 
sign decisions. It appears, nevertheless, that typical 
boundary layers may reduce the negative damping 
factor—i.e., the degree of instability ——by an order of 
magnitude (this statement applies only to long, mono- 
coque cylinders, for which the critical wavelength and 
boundary-layer thickness are of the same order of 
magnitude), thereby providing a possible explanation 
of previous discrepancies between theory and observa- 
tion. 


(2) EQUATIONS OF MOoTION—-BOUNDARY LAYER 


We require the equations of motion that govern the 
disturbed flow produced in a prescribed boundary layer 
by the surface wave displacementt 


exp { —ik [(x — Vt) cosa+ysina]} (2.1) 


of the boundary from its equilibrium position z = 29; 


t Cf. Eq. (3.13) of reference 2, x; and y; therein being replaced 


here by x and y. 
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x and y, together with z, are Cartesian coordinates in a 
reference frame fixed in the panel, ¢ is the time, V the 
phase velocity with respect to the x axis, a the angle of 
propagation (of the wave normal) with respect to the 
x axis, and k the wave number—viz., 


k = (2.2) 


A, being the wavelength of the surface wave. The 
amplitude of the wave, although for convenience posed 
as unity in Eq. (2.1), is assumed to be small compared 
with all other characteristic lengths —in particular, both 
the wavelength and the boundary-layer thickness. 
(We remark that the effect of the boundary layer on 
panel flutter almost certainly would be unimportant if 
the amplitude were of the same order as, or large com- 
pared with, the boundary-layer thickness.) 

Our model of the boundary layer is an inviscid, non- 
heat-conducting, parallel, shear flow that, in the ab- 
sence of the traveling wave disturbance, is character- 
ized by the velocity profile U/(z), the temperature pro- 
file 7(z), and the free-stream pressure ~; (which is con- 
stant through the boundary layer). The sonic velocity 
a(z) and density p(z) then are related to the temperature 
according to 


a*(z)/a> = pi/p(z) = T(s)/Ti (2.3) 


by virtue of the constant pressure; the subscript 1 on 
the fluid variables denotes the free-stream conditions 
outside the boundary layer. We neglect the x and y 
variations of all of these quantities with the implicit 
assumption that local values may be used (as in the 
boundary-layer stability problem). 

If the boundary layer is laminar, the essential ap- 
proximation in our shear flow model is that the viscous 
forces associated with the traveling wave are negligible 
compared with the corresponding inertia forces (the 
neglect of heat conduction is almost certainly of less 
importance than the neglect of viscosity). This ap- 
proximation is similar to that adopted in the inviscid 
approach to the boundary-layer stability problem,’ 
although in the present problem the approximation may 
remain valid in the neighborhood of the boundary by 
virtue of the finite acceleration there. As in the 
boundary-layer problem, the approximation breaks 
down in an inner viscous layer (see reference 3, p. 136), 
where the inertia forces tend to zero. The result is a 
singularity in the inviscid equations of motion [see 
Sections (4) and (5)], which may be said to represent 
the inner viscous layer in the limit of zero viscosity. 
We infer from the known results for the boundary-layer 
stability problem that this representation is significant 
for small but finite values of the viscosity if the dis- 
turbed motion is unstable. 

If the boundary layer is turbulent and we interpret 
U, T, p, a, and p; as mean values, the approximations 
implied by the model are more severe. First, we 
neglect the perturbation Reynolds stresses associated 
with first-order} coupling between the perturbation 


7 Cf. Appendix, reference 7. 
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flow and fluctuations in the original flow. Second, we 
neglect cross-correlations between the fluctuations of 
unlike quantities—e.g., density and velocity. We 
suggest that the latter approximation is not likely to 
introduce uncertainties greater than those inherent in 
presently available profile data for compressible bound- 
ary layers; this conjecture is supported by the observa- 
tion that these cross-correlations would not appear in 
the corresponding incompressible problem and probably 
represent only minor compressibility effects compared 
with the variation of mean Mach Number (which /s 
incorporated in our model). 

We next remark that perturbations having the phase 
velocity V’ with respect to the steady shear flow will 
appear as steady disturbances to an observer moving 
with this same velocity |’, with respect to whom the 
original flow will have the velocity distribution U(z) — 
V. This allows us to use Ward's result‘ that the per- 
turbation pressure p — /; and the perturbation velocity 
components v and w (along y and gz, respectively) may 
be derived from a potential ¢ according to 


(2.4a) 
o(U — Vv = 4, (2.4b) 
p(U — V)w = (2.4c) 
where ¢ satisfies 
(1 — M*)brr + by + O22 — 
(2/M) (M,o, + 4.9.) = (2.5) 


with J/ as the apparent (to an observer moving with 
the wave at elevation z) Mach Number 


M = M(z) = [U(z) — V]/a(z) (2.6) 
We emphasize that ¢ is not a velocity potential in the 
classical sense and that the flow is not irrotational. 

The linearized boundary condition to be imposed 


on ¢ at the wall (z = 2») is, from Eqs. (2.1) and (2.4¢e),1 


¢. = p(U — V) (Dé/Dt) = 


cosa p(U — V)%, (2.7) 


This boundary condition (together with the anticipated 
form of the expression for the pressure at z = 2) sug- 
gests that we assume a solution to Eq. (2.5) in the form 


t We note that if 1/(z) is introduced in place of 2 as an inde- 
pendent variable, the boundary condition (2.7) may be regarded 
as imposed along the streamline on which lJ = J/) = M(z») in 
the absence of the disturbance, after which the variable z may be 
interpreted as a parametric streamline coordinate, rather than 
the linear distance from the wall. This strategem avoids the 
restriction (that would be implied if s were interpreted directly 
as a linear measure) to amplitudes that are small compared with 
any distance in the boundary layer over which appreciable 
changes in velocity occur; such a restriction would be highly unde- 
sirable because of the large velocity gradient near the wall. 
We also note that, while UV = 0 at the wall, it might be desirable 


in some applications to impose the boundary condition (2.7) at 
some other streamline—e.g., the edge of the laminar sublayer. 
This would be permissible insofar as the distance between this 
streamline and the wall were negligible compared with the wave- 
length. 
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— Vt, y, 3) = pray” sec a — Vt,y) (2.8) 
Substituting this in Eqs. (2.5) and (2.7), and introducing 
m = M(z) cosa = [U(z) — V]/a(z) seca (2.9) 
we obtain 
m*(d dz) |m~*(df/dz)| — k?(1 — m?)f = (2.10) 


df 
and = —km?(z9) (2.11) 


The differential equation satisfied by f(z), to which 
function the perturbation pressure is proportional, was 
obtained by Lighthill® in his studies on shock-wave- 
boundary-layer interaction; it is appreciably simpler 
than the differential equation satisfied by w (the latter 
usually has been used in boundary-layer stability 
studies’). We also note that if m = 0(U = V) at, 
say, 3 = 2, in the boundary layer, the differential equa- 
tion (2.10) has a regular singularity there.* The im- 
mediate (zs) neighborhood of this singularity consti- 
tutes the inner viscous layer, where the inertial forces, 
being proportional to (U7 — JV)’, cannot dominate the 
viscous forces; however, the thickness of this layer ap- 
proaches zero for large Reynolds Numbers (see refer- 
ence 3, p. 136) if the disturbance is self-excited, and we 
therefore neglect it in the following development (except 
insofar as it is represented by the singularity at z = 2,). 

The differential equation satisfied by /(z) outside of 
the boundary layer, where m assumes the constant 
value m7, reduces to 


(d*f/dz*) — Bf = 0 (2.42) 


where 


8= k(1 m,2)'? = 

kil — ((U, — V)/a sec (2.13)t 
We then have f(z) ~ e (2.14) 
where z = 2, denotes the edge of the boundary layer; 
from this we deduce the boundary condition 


(df/dz) + Bf =0, (2.15) 


The finiteness condition that f be bounded as zs > © 
requires R] 8 > O, while the radiation condition re- 
quires Im £ to have the same sign as 1 — U if RIB = 0. 
These requirements will be satisfied everywhere in the 
complex |” plane if we draw the cuts for 6 from the 
branch points 1° = U, + a, sec ato + @© and exclude 
points on their top sides. Of course, it is not necessary 
to satisfy the finiteness and radiation conditions every- 
where in the I” plane, but only at those points for which 
solutions to Eq. (2.12) are to be considered (thus, in 
reference 1, the cuts were drawn into Im Il’ > O and 
only points in Im I’ < O considered). In the final 
analysis, we may deform the cuts in any way that en- 

* This singularity is only apparent if (U’/T)’ vanishes at 
s = z,; cf. Eq. (4.7). 

+ In Section (5), 8 is redefined as «(1 — m2)!/?, where x = 6, 
and 6 denotes boundary-layer thickness. 


sures the existence of the solution to the differential 
equation and the satisfaction of the boundary condi- 
tions. We remark, in this connection, that a periodic 
disturbance of the type (2.1) cannot exist by itself, 
and in a complete initial value problem we would be 
led to consider an integral superposition of such dis- 
turbances over a prescribed path in the |’ plane. 

The pressure on the panel at z = %, as obtained from 
Eqs. (2.4) and (2.8), is given by 


Pa = Pi) z=2 
— V/U1) — Vt,y) (2.16) 


where we have introduced 
F(U,/A, V/U,) = f(2) (2.17) 
and A=a, seca (2.18) 


Our boundary value problem now may be posed: 
find a solution to the differential equation (2.10), sub- 
ject to the boundary conditions (2.11) and (2.15), and 
evaluate it at z = 2. 


(3) EQUATIONS OF MOTION—PANEL 


We present in this section a brief derivation of the 
eigenvalue equation for the panel. A more detailed 
derivation has been given in reference 2. 

We assume the existence of an operator L, such that 
Lé is the reaction force per unit area associated with a 
transverse displacement ¢. L may imply both differ- 
ential and integral operations with respect to the space 
coordinates, but not with respect to the time—.e., 
Lé is in phase with ¢. The equation governing the 
panel motion then reads 


Le + /Ot?) = pi — pu (3.1) 


where /; denotes the perturbation pressure of the fluid 
beneath (inside) the panel, p, the perturbation pressure 
outside the panel, and o the panel mass per unit area. 
If the compressibility of the internal fluid is neglected, 
p; represents a virtual inertia force, which, for the 
traveling wave defined by Eqs. (2.1) and (2.2), is given 
by 


Pi = pik(V cos (3.2) 


where p; is the fluid density. (This follows directly 
from the fact that the effective depth of an incompres- 
sible fluid for a surface wave motion is | /k; see reference 
2 for an alternative derivation.) 

The eigenvalue equation corresponding to the travel- 
ing wave of Eq. (2.1) may be derived (see reference 2 
for details) by relating the operator L to the wave 
speed lV, that would result in the absence of the aero- 
dynamic force (but including the virtual inertia force 
of the internal fluid); thus, substituting ¢ and p, from 
Eqs. (2.1) and (3.2) in Eq. (3.1) with V = V% therein, 
we obtain 


Le = [o + (pi/k)] cos (3.3) 


We emphasize that lV’) is an eigenvalue of Eq. (3.3) 
and is a function of k and a. 


ve 
of 4 
ve 
to 
in 
a- 
in 
ly 
| 
se 
Ig 
1e | 
y 
y 
) 
)) 
)) 
h 
) 
ig 
1 
1 
4 


84 


Substituting p, from Eq. (2.16), and L¢ from Eq. 
(3.3) in the complete equation (3.1), with V # Vj, in¢ 
and p;, and dividing through by the coefficient of V, 
we obtain the eigenvalue equation (for V), 
V2? + pA2F(U,/A, V/Ui) = Vo? 

pi/(pi + ok) 


is the aerodynamic-mass parameter, and F and 4 are 
defined by Eqs. (2.17) and (2.18). 


(3.4) 


where (3.5) 


(4) SomE GENERAL CONSIDERATIONS REGARDING 
STABILITY 


The explicit determination of the stability boundaries 
for a given panel, as represented in the specification 
or derivation from Eq. (3.3)] of Vo as a function of k 
and a and the mass parameter yp, and a given boundary 
layer, as represented by the specification of velocity and 
temperature profiles, requires the simultaneous solution 
of the differential equation (2.10), subject to the bound- 
ary conditions (2.11) and (2.15), and the transcendental 
equation (3.4).* This represents a formidable prob- 
lem, which is especially complicated if V = U within 
the boundary layer, giving rise to a singularity of the 
differential equation (2.10). In order to clarify the 
role of this singularity with respect to panel flutter, 
we develop in this section an implicit solution to the 
boundary value problem. 

We first rewrite Eq. (3.4) in the form 


V = Voll — w(A/Vo)? F(Ui/A, V/U)}"? (4.1) 


where the phase of the radical is defined such that 1) = 
Vo at w = 0. It then follows, neglecting structural 
damping (Vo real), that the imaginary part of |’ will 
be opposite in sign to the imaginary part of F, so that 
instability requires Im F > 0. 

We obtain an implicit solution for F by multiplying 
Eq. (2.10) through by m~*f (where the bar denotes the 
complex conjugate), integrating from z = 2 toz = 2, 
integrating the term f(m~*f’)’ by parts, imposing the 
boundary conditions (2.11) and (2.15) to eliminate f’ 
at the limits, and solving for f(z); the end result is 


F = f(z) = (1 — m2)? + 
[ 24 (1 m?) If 


We shall assume that the phase angle of V is small (V 
approximately real); then m will be approximately 
real, and the imaginary parts of the first and second 
terms on the right-hand side of Eq. (4.2) must be de- 
rived chiefly from the radical (1 — m,?)'/? and the singu- 
larity at m = O, respectively. 

The phase of (1 — m,2)'”? for V approximately real 
is, by definition, 0 for im, | < lor +7/2 for +m, > 1, 
corresponding to our previous choice of branch cuts 


2] m—*dz (4.2) 


* Alternatively, Eq. (3.4) may be combined with Eq. (2.11) 
to obtain a homogeneous boundary condition on f at s = 2. 
This formation leads to a more symmetric boundary value prob- 
lem but is less satisfactory for the approximate solutions treated 


herein. 
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for [see remarks following Eq. (2.15) |. Then, writing 


F, = f(a) [2 (m2 — 1)", 


we infer that the contribution to Im F(= —Im F) of 
the first term in Eq. (4.2) is F,, 0, or — F, for V < U; — 
A,U,-A<V<U,+A,o0r V> U, + A, respec- 
tively, if V lies close to the real axis but not too close 


| >1 (4.8) 


to U; A. The integral in Eq. (4.2) will contribute 
an imaginary part, say F,, to F only if 0 < V < U,, 
in which case m will vanish at the point z = 2. It 
follows that 
Im 0<V<U,-A 

= F,, U%,-A<V<U, 

0, 


- F, U,+ A <7 


We remark that if V does approximate U, + A, as in 
Section (6), the contribution of the first term in Eq. 
(4.2) to Im F must be modified, but this does not 
affect the contribution of the second term, which re- 
mains as F, or 0. 

We may deduce a simplified expression for /, by 
not ng from the differential equation (2.10) that f’(z) = 
O(m) and f(z) = f(z.) + O(m*) as zs > 2,;f accordingly, 


F, = —k-'Im |? + 
k? (1 — m’) if m-*dz (4.5a) 
= —k6| f(z.) |? Im K; (4.5b) 
where K,, as defined by 
K, = 6“! { (m-* — 1)dz (4.6) 


is the integral appearing in the solution for /\(z) dis- 
cussed in Appendix (A). The sign of Im A, evidently 
depends on how the singularity is circumvented; this is a 
delicate question, but it has been examined thoroughly 
by Lin (reference 3, ch. 8), and we merely state his 
result that the path of integration must be indented 
over the point z = z,,{ whence 
Fy = X 

(d/dz)(U'/T) (4.7) 


[see 


— tka,” sec? a | (Ze) 


where the subscript c implies evaluation at z = 
Appendix (B) for details]. 

The integral A, occurred in the laminar boundary- 
layer stability problem, and Lees and Lin*® have shown 
that if the original flow absorbs energy from a small 
disturbance in the boundary-layer profile 


(d/dz) (U'/T) | < 0 (4.8) 
It does not follow that Eq. (4.8) is a necessary and suffi- 


cient condition for the stability of the boundary layer, 


Cf. Eq. (5.4b). 

t We note that the sign of the wave number usually assumed 
in boundary-layer studies is opposite that chosen here, which has 
the effect of interchanging 7 and —7; thus, the path of integration 
elected by Lin actually is indented under the singular point 
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but it seems likely [on the basis of empirical consider- 
ations, as well as Eq. (4.8)] that it will be satisfied for 
both laminar and turbulent boundary layers in those 
aerodynamic regimes for which panel flutter presents 
a serious problem. We conclude, therefore, that F, 
will be positive and destabilizing. 

A rough approximation to the value of z, for a turbu- 
lent boundary layer may be deduced from the well- 
known one-seventh power law—viz., 


= [(z — 2%)/6]'”” (4.9) 


where 6 is the boundary-layer thickness. We antici- 
pate that panel flutter will be most serious when V 
approximates or is less than U/; — a [as in the absence 
of a boundary layer; ef. reference 2 and Sections (6) and 


(7) below]. Setting U = U, — a and z = ¢, in Eq. 
(4.9) then yields 
— 2)/6 = [1 — (4.10) 


which has the values 5 X 10-4 and 8 X 10~* for MW, = 
1.5 and 2, respectively; these values are so small that 
they might place zg, in the laminar sublayer, where the 
one-seventh power law no longer would be valid, but 
they may be assumed to give the correct order of magni- 
tude. Of course, (3, — 2)/6, as given by Eq. (4.10), 
does not remain small for large 1/,, but it seems likely 
that F, then would be unimportant compared with F,. 
We add that (zs, — 2)/6 probably would approximate 
(1 — 1 J) for a laminar boundary layer, but a laminar 
layer almost certainly would be so thin as to have no 
effect on panel flutter. 

The negative damping represented by F, is associated 
with coupling between the surface wave on the panel 
and a sound wave in the air; it is qualitatively similar 
to the negative damping predicted in the absence of the 
boundary layer and is discussed in more detail in refer- 
ence 2. The presence of the boundary layer reduces 
F, by reducing the acoustic coupling between surface 
wave and sound wave, but the energy transferred to 
the panel is derived primarily from the outer flow. 
The negative damping represented by F,, on the other 
hand, is associated with coupling between the surface 
wave and the shear flow in the boundary layer and 
could not arise in the absence of the boundary layer 
(more precisely, F,, tends to zero with ké). The cor- 
responding energy transferred to the panel is derived 
entirely from the shear flow through a mechanism 
basically similar to that found in (the inviscid approxi- 
mation to) classical examples of hydrodynamic insta- 
bility. It is possible that the boundary layer could be 
responsible for panel flutter if — A < Vinin < 
but if Vinin << Uy; — A, its net effect is more likely to be 
stabilizing in consequence of the reduction of F,,* [but 
see the penultimate paragraph of Section (6) }. 


* We note that F, = 0 in subsonic flow, and instability then 
would have to be charged entirely to the shear flow in the bound- 
ary layer. It seems likely that F, would be much too small to 
render practical aircraft structures unstable, but it may afford an 
explanation of the generation of gravity waves on a liquid 
surface .7 
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(5) THe BouNDARY VALUE PROBLEM 


We now return to the boundary value problem posed 
by the differential equation (2.10) and the associated 
boundary conditions (2.11) and (2.15). Introducing 
the dimensionless variable 


= (s — (5.1a) 
and the boundary-layer thickness parameter 

= kb = (5. 1b) 
where 6 = — (5.1¢) 


in Eqs. (2.10), (2.11), and (2.15), we obtain 
m*(m~*f’)’ — «(1 — = 0, O< F< 1 (5.2) 
fo’ = (5.3a) 
fi’ + Bfi = 0 (5.3b) 


where the primes imply differentiation with respect to 
¢, the subscripts 0 and 1 imply evaluation at ¢ = 0 and 
1, and @ is defined by Eq. (2.13) with k now replaced 
by « therein. 

The differential equation (5.2) has a regular singu- 
larity at ¢ = ¢,, where m(f,) = 0. The exponents of 
this singularity are 3 and 0, and we find that the cor- 
responding solutions in its neighborhood are of the 
form 


= — — (3g/4) (F — + 
— (5.4a) 


= log [(¢ — + 
1 — (x?/2) (§ — + OF — (5.4b) 
where = (5.5) 


The solution w; is everywhere real for real values of 
¢ — ¢., while w. has been normalized to have the imagi- 
nary part O or (— rigx?/3)w, as < or > respec- 
tively (since the logarithmic branch point must be 
circumvented by indenting over ¢ = ¢,.). The Wron- 
skian of this pair of solutions is 


Wi wr, wet = — wow,’ = (5.6) 


We also find it convenient to introduce a second pair 
of solutions, f; and fs, that satisfy the boundary condi- 
tions 


fio = 1 (5.7a) 
fo’ = 0 (5.7b) 
foo = O (5.7¢) 
foo’ = (5.7d) 


where the first and second subscripts identify the solu- 
tion and the point of evaluation, respectively. The 
corresponding Wronskian is 


Whfy fol = (5.8) 


Invoking the boundary conditions (5.7) on linear com- 
binations of w,; and w. and making use of the Wronskian 
(5.6), we obtain 
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Fic. l(a). The real part of m,2K,, as given by Eqs. (B-16)-( B-20) 
for a = 0 and an insulated boundary layer. 


and 
= [we — Ww (5.9b) 


The required solution, subject to the boundary con- 
ditions (5.3), now may be expressed in terms of either 
w, and or f; and fz. We find 


+ Bw) Woy’ — + Who" 
(5.10a) 


and f = FACE) — Khoo) (5.10b) 
where F = x[(fa’ + + | (5.11) 


We require only the imaginary part of F throughout 
much of the subsequent analysis. This can be calcu- 
lated directly from Eq. (5.11), but it may be advanta- 
geous to consider separately the components F, and F,, 
as defined by Eqs. (4.3) and (4.5b), for the results then 
do not depend on the direct calculation of phase (which 
might be rather small in some applications). Setting 
¢ = 1 in Eq. (5.10b), simplifying the result with the 
aid of the Wronskian (5.8) evaluated at ¢ = 1, and 
substituting in Eq. (4.3), we obtain 


F, = — 1)'?/\ ful + Bhu 


Setting ¢ = ¢.in Eq. (5.10b), expressing the denomi- 
nator of the result in terms of f; with the aid of Eq. 
(5.9a), and substituting in Eq. (4.5b), we obtain 


(5:12) 
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+ Bwn |? 
9 | far’ + Bfu |? 


The solutions fi, f2, and w, are developed as power 


| Im A, (5.13) 


series in « in Appendix (A). The first approximations 
to the quantities required for the evaluation of F, and 
F, are 


fu = 1 + O(k?) (5. 14a) 
fu’ = [1 + O(x?)] (5.14b) 
wu = [L + (5. 15a) 


Wy! = (3m,2/m,"2) [1 + O(«?) ] (5.15b) 
where A, is the integral introduced in Eq. (4.6) and | 
evaluated in Appendix (B), and J; is defined by Eq. 
(A-15) in Appendix (A). Substituting these approxi- | 


mations in Eqs. (5.12) and (5.13) yields : 
— 
F, = + Of«*)}, 
|—i(m?2 — 1)? + «m)?K, 
m, > (5.16a) 
= 0 m <1 = (5.16b) 
and 
1/2 
my? (1 — Jy? 
(1 — + 


Im K, [1 + O(x?)]  (5.17a) 
Im Ky 
9\ 1/2 Xx 
(1 — + |? 
[1 + O(x?, — m,")| (5.17b) 


The approximation (5.17b) anticipates (on the basis of | 
the results for « = 0; see reference 2) that supersonic 
panel flutter should prove most serious for values of 
m, only slightly larger than unity for sufficiently small 
values of x. 

The approximations to F, and F, provided by Eqs. 
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Fic. l(b). The imaginary part of —mm,*A), as given by Eqs. 
(B-16)-(B-20) for a = 0 and an insulated boundary layer. 
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Fic. 2(a). Fa versus V/U; for M, = 1.5. 
S. 
(5.16) and (5.17b), together with the values of A, 
“ given by Figs. 1(a) and 1(b) (for a = 0 and an insulated 
| | boundary) are plotted in Figs. 2(a, b) and 3(a, b) for 
M, = 1.5 and 2.0, and «x = 0.5, 1.0, and 2.0.* Also 
4 plotted (labelled ‘‘exact’’) are some results that were 
|| obtained by numerical integration of the differential 
| equation (5.2) on a digital computer.f The agreement 
— 
| *More extensive results are plotted in Ramo-Wooldridge 
Report GM-TR-299. 
| + The procedure was to obtain first the regular solution w, by 
integrating away from the singularity, starting the solution by 
power series expansions near ¢ = ¢, and then continuing by the 
Rungé-Kutta method over the range (27, 1). The singular 
solution w. then was expressed as 
= (qx?/3)ar log — + 
- the inhomogeneous differential equation for g determined by 
substitution in the original differential equation, and the regular 
function g obtained by the procedure outlined for w;. In those 
s cases where ¢, was outside of the range of integration, f; was 
determined directly by integrating out from ¢ = 2714. The de- 
tails and programming were worked out by Dr. Samuel Conte 
qs: 


and David Bussard; the computer was a Remington-Rand, 
Univae Scientific Model 1103A. 
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doubtless could be improved by including terms of 
higher order in x in the approximate results for F, and 
F, or, in some cases, merely by using Eq. (5.17a) in 
place of Eq. (5.17b), but the errors in the approxima- 
tions (5.16a) and (5.17b) are not likely to exceed those 
already present in our simplified model of the boundary 
layer. 

We emphasize that the results (both ‘exact’? and 
“approximate’’) for F, are not valid for 1” U less than 
about 0.3 in consequence of the breakdown of the one- 
seventh power law assumed in the calculations. An 
adequate correction factor for F, could be achieved by 
forming the ratio of Im A, for the actual U’(s) and 7(z) 
to that based on the one-seventh power law, using Eq. 
(B-Sc). The results for /, do not appear to be very 
sensitive to changes in the profile. 


(6) APPROXIMATE STABILITY ANALYSIS 


We consider in this section the approximate calcula- 
tion of the negative damping ratio (in the notation of 
reference 2) 


0 0.1 0.2 #043 04 O14 0.6 0.7 0.8 0.9 1.0 


v/U, 


Fis. 2(b). Fa versus V/U, for M, = 2.0. 


5) = —Im V/RIV (6.1) 


on the assumption that the second term in the braces 
of Eq. (4.1) is sufficiently small compared with unity 
to permit the approximation 


V = Voll — (1/2)u (A/Vo)? F(U/A, Vo/Ui)} (6.2) 


Eq. (6.2) may be regarded as the first approximation 
based on a power series expansion in the parameter uy, 
and it will be valid for » < 1 insofar as F remains 
bounded. Now [see Eq. (3.5) | almost certainly would 
be small in practical applications, while |F can be 
large only if both |V — (U, + A)| < U, and x <1; 
the latter contingency is examined in the following 
section, but we remark here that it is apt to be impor- 
tant only when « is of the same order of magnitude as 

1/3 

Substituting Eq. (6.2) in Eq. (6.1) and evaluating Im 
F from Eq. (4.4) yields 


= (1/2)u(A/Vo)? (Fa + Fy), O<Vo<U,-—A (6.3) 


in the region of principal interest. The maximum value 
of 69 with respect to V> almost certainly will occur in 
this region if (Vo)min < Ui — A, but F, = Oif U — 
A < Vo < Ui; if Vo > U, panel flutter of the type con- 
sidered here is not possible. 

Perhaps the most expedient approach to the deter- 
mination of the maximum value of 6) with respect to 
variations of x and a@ is to adopt the approximations 
of Eqs. (5.16) and (5.17b) for F, and F,, determine the 
corresponding value of (6))mar, and then estimate the 


—— Equatione (5.16) and (5.17b) 


IN core exact 
0.24 
0.20 | 
\ 
% 
Fy \ 
8:0 
0.12 N 
\ 
“J \ = 2.0 
0.08 
/ 
\ 
\ 
0.04 K 
\ 


Fic. 3(a). F, versus V/U; for M, = 1.5. 
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approximate corrections to F,, and F, from the results 
of Fig. 2. Substituting Eqs. (5.16) and (5.17b) in 
Eq. (6.3) yields 


1 A \? 
bo = m,” Im (6.4a) 
2 Vo — + «m2 


— 1 — xm,’ Im Kj] 
(Vim? — 1 — Im + [xm,? RI K,]}? 
m, > (6.4b) 


We recall that [Eq. (2.9) | 


m, = [(U, — Vo)/a:] cosa = 
M, [1 — cos (6.5) 


while A, depends on Vo, Ui, /;, and cos a; as suggested 
in Appendix (B), the error in assuming that A, depends 
only on 1, cos a, rather than .1/; and cos a separately, 
usually will be small (the essential approximation is | 
the replacement of 1/, by 4; cos a in the temperature 
profile). 

We consider as a more specific example a pressurized © 
cylindrical shell for which the internal pressure is sufli- | 
ciently high to render negligible the effects of bending 
on J, which then is given by [see reference 2, Eq. | 
(5.5), with D = Oand N, = 2N, = p:R| 


= (1/a) { Eh(cos a/kR)? + 
+ tan? (6.6) 


—— Equations (5.16) 
and (5.17) 


exact 


Fic. 3(b). F, versus V/U, for MM, = 2.0. 
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We introduce as the independent variables 


x= (6.7a) 
y = cos? a (6.7b) 
and, as a measure of the wavelength, 
z= 2,/kR (6.8a) 
= V2Eh/p,R (6.8b) 


Substituting Eqs. (6.7a, b) and (6.8a, b) in Eq. (6.6), 
we obtain 


2 = V(1/y) [1 + (x/xm)?] — (2/92) (6.9) 


where Xn = WV p,R/20U? (6.10) 
denotes the minimum value of x (v = x, at a = Oand 
k = ©). The boundary-layer thickness and aero- 


dynamic mass parameters [the latter given by Eq. 


(3.5) with p; = O therein] then may be expressed as 


k = = 2, (6/R)2— (6.11) 

and pi/ok = (6.12) 
Noting also that 

= — x)*y (6.13) 


and, on the basis of the aforementioned assumption 
for Ki, 


= Ui, M cos a) 


(6.14) 
Eq. (6.4b) may be transformed to 
bo = ( 2)(pR o2)2[(1 x 
Im! —iV — x)*y? — 1+ 
(x, My'”)}= (6.15) 


We require the maximum value of 60, as given by Eq. 
(6.15) in conjunction with Eq. (6.9), with respect to 
independent variations of x and y over the intervals 


Ss 4 (6.14a) 
2/[1 + (x/xm)?] Sy 1 (6.14b) 


The required calculations would be difficult to carry 
out analytically, but they are quite simple on a high- 
speed computer. 

If the boundary layer is relatively thin (x < 1) we 
find that (6)maz is likely to occur for values of x and y 
such that mz, exceeds unity by only a small amount and 
the two terms in the braces of Eq. (6.15) are approxi- 
mately equal in magnitude. The imaginary part of 
K,—and, therefore, the negative damping effect of 
profile curvature—then is likely to prove rather unim- 
portant, although its real part plays an essential role 
in preventing the infinity that would be indicated by 
Eq. (6.15) at m, = 1 if « — O* (in which case the 
analysis of the following section would be required). 

* The effect of the boundary layer then may be compared 
with that of small damping in limiting the response of a simple 
oscillator at resonance. 


The approximations of Eqs. (5.16) and (5.17b) evi- 
dently are entirely adequate for this case. 

If the boundary layer is relatively thick (say « > 1 
or 2) and x,, is not too small, it appears that (60) mer 
is likely to occur at m, = landa=0. The instability 
in this case would be associated entirely with profile 
curvature, and it therefore would be important to im- 
prove the approximation to F,, relative to that of Eq. 
(5.17b), and to investigate the effects of possible 
departures from the one-seventh power law. 

Calculations based on the results of this section have 
been carried out for pressurized, monocoque shells 
filled with either gas or liquid and compared with cal- 
culations based on the results of reference 2. It ap- 
pears that the degree of instability for supersonic flutter 
may be reduced by an order of magnitude (roughly a 
factor of 5-10) in consequence of boundary layers for 
which x > 1/2. 


(7) STABILITY ANALYSIS FOR A VERY THIN BOUNDARY 
LAYER 


The analysis of the preceding section breaks down 
in the neighborhood of m, = 1 if x is of the same order 
as u'’* (of course, if « is sufficiently small compared 
with yw the effects of the boundary layer on panel 
flutter will be negligible, but turbulent boundary layers 
are not apt to be this thin at supersonic speeds). The 
assumption that both « and im — 1 | are small allows 
us to approximate F by [see Eqs. (5.11) and (A-Sa)- 
(A-9b); also ef. Eqs. (5.16) and (5.17a), with which this 
approximation is consistent | 


F = — m)'? + X 
(1 + O(«2)] (7.1) 


Substituting F (from (7.1) in (3.4) and eliminating mm, 
through (2.9) and (2.18), we obtain, after some alge- 
braic manipulation, 


[(V? — Vo?)/A2] — — V)/AP}'? + 
— V)/A]? V/U,)) + 
— V)/A}? =0 (7.2 


where the dependence of A, on both l’/4 and VV’, U; has 
been explicitly denoted. 

We next introduce the dimensionless variables v and 
representing the departures of and from U, — 
A, according to (cf. reference 2) 


V = (U, — A) [1 + (1/2) (M, — 
—3n/2<argy < 7/2 (7.3a) 

Vo = (Ui — A) [1 + 
(7.3b) 


where the permissible range of arg v is deduced from the 


restrictions on V and 8 [Eq. (2.13) ], is real (structural 
damping being neglected), and .\/, denotes the (free- 
stream) Mach Number relative to an observer moving 
with the wave front—viz., 


M, = U,/A = U, cos (7.4) 


in 
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Fic. 4. h(n), the ratio of the maximum negative damping ratio 


with boundary layer to that without. 
We also find it convenient to introduce the boundary- 
layer parameters 6 and y according to 


@ = («/3) {2[(M. — 1)/u]}'*K, 1 — 
(7.5a) 
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n = 3'6 (7.5b) 


We will use the parameter @ initially, but we find » 
more convenient in the end results. Substituting Eqs. 
(7.3a, b) and (7.5a) in Eq. (7.2) and neglecting higher 
powers of u then yields 


(v — 3&) (v'/? + 30) +2 =0 


We observe that the importance of the boundary layer 

in this approximation is determined by the size of 3 6 

relative to v'’*; in particular, «x < 1 is not a sufficient 

condition for the neglect of the boundary layer. 
Rewriting Eq. (7.6) in the form 


(7.6) 


(v'/? + 6)3 — + 6) + (A* + B®) = 0 (7.7) 
where 
= {(1 — + 6%) + 

(1 — + 6%)? — + (7.8) 
wehave 4e+ & = (7.9) 


where the three roots of the cubic equation (7.7) corre- 
spond to the three cube roots (€) of —1 (é€ is the com- 
plex conjugate of €). This result is valid for complex 
values of both é and 6, but the explicit determination of 
the imaginary part of v, on which the damping ratio 
depends (assuming » to be approximately real), is pos- 
sible only if both and @ are assumed to be real; on this 
last assumption we find 


Im = —(9/3/2) | — + B— 201, 
(1 — 3&6 + 6%) > + 6%)? (7.10a) 
= 0, (1 — 38 + 6) < + 6)? (7.10b) 


If 6 is not real Im vy must be determined numerically 
from Eq. (7.9). 

The maximum value of —Im v with respect to the 
spectral parameter & is determined in Appendix (C). 
We designate the ratio of this maximum to its value 
in the absence at the boundary layer (3! 2/2' ata =0 
and — = 0) as h(n), which is given by Eqs. (C-7) and 
(C-S) and is plotted in Fig. 4; &,,, the value of £ at which 
this maximum occurs, is plotted in Fig. 5. The maxi- 
mum (with respect to &) negative damping ratio then 
is given by [ef. Eq. (4.13), reference 2| 


5) = (—Im V/RIV);-;,, = 


— (7.11) 


We emphasize that the results for 7 < 0 may not be 
directly significant, since A; and 7 will be positive real 
if the singularity at z = 2, is excluded, while if it is in- 
cluded A, and » must be complex; however, these re- 
sults may be useful in determining a first approximation 
to 6) for complex 7. 

We find, by numerical comparison, that the approxi- 
mations £, = —7? and h = 3n) are adequate for 
n > 1 [Eqs. (C-10) and (C-11)]. Substituting these 


approximations in Eqs. (7.3b) and (7.11) and elimi- 
nating » via Eq. (7.5) yields 
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5) = wl4eK,(M, — 1)2]-, (7.12) 


Volemtm = (Ui — A) — (1/2) («K;)*A, 


9 > 1 


n>1 (7.13) 


The real and imaginary parts of A,, as evaluated in 
Appendix (B), are plotted in Fig. 6 for V/U; = 1 — 
M,~'! (ef. Fig. 1). As suggested above, the error in 
5) | replacing WW, by MW, cos a for a ¥ 0 is not likely to be 


er large. 
(S) CONCLUSIONS 


We conclude that the presence of a boundary layer 
has two distinct effects on the flutter of an infinite 
panel (infinite being an aerodynamic approximation 
implying wavelengths small compared with the dis- 
tances between constraining supports, as for a mono- 
coque shell). First, the acoustic coupling between 
surface waves on the panel and sound waves in the 
outer flow (this coupling can transfer energy to the 
panel if the surface wave speed is less than the difference 
9) — between flight speed and sonic speed) is reduced. 
Second, energy can be transferred from the shear flow 
“| in the boundary layer to the panel if the surface wave 
speed is equal to the local speed in the shear profile 
at a point of strong profile curvature. Our simplified 
analysis (neglecting the direct effects of viscous forces 
and turbulent fluctuations) indicates that the former 
effect should predominate in typical configurations 
and may reduce by an order of magnitude the degree of 
instability predicted in the absence of the boundary 
layer. We emphasize that our results should be re- 
garded primarily as contributing to the basic under- 
standing of panel flutter; but, in the absence of experi- 
) | mental confirmation, they should not be regarded as 
| quantitatively adequate for important design decisions. 


APPENDIX (A)—SOLUTION OF DIFFERENTIAL EQUATION 


We require solutions to Eqs. (5.7) and (5.8)—viz., 


m*(d/d¢) [m-*(df/dg)| — — m*)f = 0, 

§ 0<¢ <1 (A-1) 
fo =1, fu’ = 0 (A-2a) 
foo = 0, foo’ = mi? (A-2b) 


_ where the second subscript indicates the value of ¢. 
A solution that is useful for small or moderate values 


may be obtained by introducing the expansions 
n=0 
- — in Eq. (A-1), which then yields 
(d/dg) [m-? (df /dg)] = (m-? — (A-4) 
We may satisfy the boundary conditions of Eq. (A-2) 
with the zeroth approximations 
fy = (A-5a) 
0 
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The succeeding terms (m > 1) in Eq. (A-3) then must 


satisfy the homogeneous boundary conditions 


ff (0) fi” (0) (A-6) 


Integrating Eq. (A-4) subject to these conditions, we 
obtain the recursion formula 


= f (m-? — 1)f"~” dé, 
0 0 


n> 1 (A-7) 


We require, in the end [ef. Eq. (2.24) ], the values of 
ti, fi’, fo, and f’ only at ¢ = 1. Following Lees and 
Lin®* (but with slightly revised notation), we may 
express these in a series of definite integrals according to 


fr = (A-Sa) 
n=0 
fu’ = > (A-S8b) 
n=0 
+4 
A 
-ImK, 
+2 
/ +RIK, 
+1 
K, 0 
-1 
\ 
-3 
-4 
2 3 + 
Fic. 6. The real and imaginary parts of K, for the special case 
a = 0(V = U,; — a;) and an insulated wall with y = 1.4. 


= 
€ 
| 


Sor’ = > (A-9a) 
n=0 
= Toy (A-9b) 
n=0 
Ko = 1 
1 
Ki = f (m-? — 
0 
1 
Ky, = f (m~* — 1)d¢é 
0 0 
1 
K; (m-2 — dar mde X 
0 0 
(m-? — 1)dg,... (A-10) 
0 
= 1 
1 
H, = 
0 
Hy, = fom nar mat,... 
0 0 


The numerical evaluation of these integrals for a /aminar 
boundary layer has been discussed in some detail by 
Lees.’ Appendix (B) contains a discussion of the 
evaluation of A, (the most important of the integrals) 
for a turbulent boundary layer. 

The analytic solution w,, as introduced in Section 
(5), also may be expanded according to Eq. (A-3). 
Integrating Eq. (A-4) subject to the boundary condi- 


tions 


we obtain w,“ = (3/m,’?) (A-13a) 


= mde (m-* — (A-13b) 
fe 


We require w,(1) and w,’(1), which may be placed in 
the forms 


wy’ = (3/m,'2) Son 
0 
wu = (A-14b) 
0 
Jo ) 
1 
J; = 
Jo = (m-? — A md¢ (A-15) 
J3 = (m-? — X 
Fe Fe 


mrdt,... | 
te ) 
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We consider the integral 
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K, = sf (m-? — 1)dz, m/(z) = 
0 
[U(s) — V]/a(z) see a (B-1) 
Introducing the dimensionless notation 


= a*(z)/a,?_ (B-2) 


where 7(é) is the relative temperature distribution, ex- 
pressed as a function of the relative velocity, A, be- 


comes 
= —1 + sec’ a (B-3) 
where 
1 
I(é.) = [rdf /(— — &,)?] (B-4a) 
0 
1 
J [r(&) ¢’(é)dE/(E — &-)?] (B-4b) 
0 
and we have assumed U = U\(é = 1) at z = 6(¢ = 1). 


Integrating Eq. (B-4b) once by parts, along a path 
indented over the singularity at ¢ = ¢,, and then sep- 
arating out the singular portion (/,), we obtain 


IT = — (1 — &) + 
B-5) 


where 


1 
(757) 0" [\dé/(—E — (B-6) 
0 


and the subscript c implies evaluation at the singular 
point & = &. Integrating /, over the singular point 
yields directly 

I, = (rt’),' fe | (B-7) 


Substituting this result in Eq. (B-1) yields for the 
imaginary part of A, 
Im A, = —7M,~ sec? a(7¢’),’ (B-Sa) 

(d/dz) (B-Sb) 


w5—' a,? sec? a[T(z,)/U" (z.)] X 


(d/dz) [U'(z)/T(s) (B-8e) 


We consider, as a particular example, the turbulent © 


velocity profile 
f=? (B-9) 
together with the quadratic (in velocity) temperature 
profile® 
r(é) =at+ + a = T/T, 
c= — 1)/2]) M2, b=1-a- cf 


for which = 7(6a& + + (B-11) 


We note that Eq. (B-9) is not accurate right down to 
¢ = 0, where it gives ¢’(0) = 0, but the error in equat- 
ing the first term in Eq. (B-5) to zero for large free- 


~ 


whet 


Subs 
from 


kK, = 


Th 
culat 


mK 


whicl 


mK 


with 


The 1 
Figs. 
insula 
with 
with « 


of MW, 


We 


ing the 


It ise 
spond; 

* The 
GM-TEI 


92 
= 
stre: 
2 x] 
neg 
in E 
= 
[= 
x 
q 
| 
_ 
| 
with r 
| 


-1) 


)) 


PANEL FLUTTER IN THE PRESENCE OF A BOUNDARY LAYER 93 


stream Reynolds Numbers (say greater than 10°) is 
negligible. Substituting Eqs. (B-9) through (B-11) 
in Eq. (B-5) and carrying out the integrations yields 


I = Tlals(é.) + ble(&) + — (1 — (B-12) 


where 
n—1 m 


1,(é-) = (n + 1) 


m=on —™m 


Substituting Eq. (B-12) in Eq. (B-3) and a, b, and c 


from Eq. (B-10) then yields 


K, = —|]+ iM, “2 sec? a} 
(1 — &)—"] + (70/71) — + 


The quantity actually required in the stability cal- 
culations is 


mK, = [(U; — Vy) Ki = 
M2? cos? a: (1 — &)?Ay (B-15) 
which we write in the form 
m?Ky = —m? + (10/11) + + 
[(y — 1)/2]M? C(é.) 
with 
= 7(1 — &)? — (B-17)* 
B(é.) = 7(1 — &) — &) — 1] (B-18)* 
C(é-) = 7(1 — — Ix(E-)] (B-19)* 


The real and imaginary parts of m,°A, are plotted in 
Figs. 1(a) and 1(b) fora = 0, AJ; = 1.5 and 2.0, and an 
insulated boundary—1i.e., one for which 


T)/T, = 1+ — 1) /2] = 1+ (1/5) (B-20) 


with y = 1.4. The latter results also may be used, 
with only a small error, for a ¥ 0 if the nominal value 
of 7, on the curves is replaced by 1, cos a. 


APPENDIX (C)—MaAxImMuM VALUE OF Im vy 


We wish to maximize the quantity 
—Im v = (v/3/2)|A — B|-|A + B — 20| (C-1) 


with respect to &, A and B being given by 


= {(1 — 36+ 6) + 
[(1 — 366 + 6%)? — + (C-2) 
Differentiating Eq. (C-1) with respect to £ and equat- 
ing the result to zero, we obtain 


+ B)@ + AB] (A + B — 26) = 0 (C-3) 


It is evident that the maximum value of —Im » corre- 
sponds to the vanishing of the first factor in Eq. (C-3); 


* These functions are plotted in Ramo-Wooldridge report 
GM-TR-299. 


dividing this through by 4Bé@ and combining the re- 
sulting equation with the values of 1B and 4’ + B*, 
obtained from Eq. (C-2), we have 
(C-4a) 
AB — (€+ &) = 0 (C-4b) 
+ — 2(1 — 3% =0 (C-4e) 
Eliminating A and B from Eqs. (C-4a, b, c) yields the 
cubic equation 
— + 29% = 0 (C-5) 
for the required value of &, say &,,. The corresponding 
value of —Im v may be placed in the form 
(—Im v) maz = = (C-6) 
where /, normalized to have the value unity at 7 = 0, 
is given by 
h = + 9? 
— 9°) [E+ (C-7) 
and 7 is introduced in place of @ because of the manner 
in which it enters & (see below). 

Eq. (C-5) has only one real root, opposite in sign to 
6, if |n| = |3'*0| < 1, while if In| > 1 there are three 
real roots; solving for these roots yields 

[1 — (1 — In| (C-8a) 


= + 0 , >1  (C-Sb) 
—2r/3 
where cos = | (C-9) 


If ‘n| > 1 that value of £ yielding the maximum value 
of his to be chosen. 
The limiting values of £ for large |n| are 


Em = — (C-10) 
The corresponding asymptotic forms of / are 
(C-1la) 
(C-11b) 


h~ 2"? /3n, 1 


We note that these approximations are quite adequate 
for |n| > 1. 
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Dynamic Longitudinal Stability Equations for 


the Re-Entry Ballistic Missile 


E. V. LAITONE* 
University of California at Berkeley 


SUMMARY 


The limitations of the linearized equations predicting the dy- 
namic stability during rapid deceleration of a nonrolling missile 
are discussed for the usual ballistic missile which has a longi- 
tudinal plane of symmetry and is trimmed to fly at a nearly zero 
angle of attack trajectory. 

In order to demonstrate that only the first-order terms can be 
rigorously justified, the ordinary second-order linear differential 
equation predicting the oscillations in angle of attack is derived 
by two entirely different approaches. The first method follows 
the usual procedure of small perturbations in the flight trajectory 
equations, while the second method utilizes Euler’s dynamic equa- 
tions, for axes that are rigidly fixed in the moving body. 

It is found that acceleration increases the damping while de- 
celeration decreases the damping of any stable oscillation. 


SYMBOLS 


(u, w) = linear velocity components of translational 
velocity vector V, along (x, z) axis, re- 
spectively (see Fig. 2) 

q = d0/dt = angular velocity component about y axis 

f = mass moment of inertia about y axis at 
center of gravity = S (x? + 2*)\dm 

(AS2) = forces along (x, z) axis, respectively (see 
Fig. 2) 

M = aerodynamic moment about y axis at center 
of gravity 

m = W/g = mass of body 

0 = angle of rotation about y axis 

a= w/u = angle of attack (see Fig. 2) 


oY = flight path angle relative to horizontal 
mass density of atmosphere 


p = 

uw = pSL/2m = nondimensional ratio proportional to mass 
of displaced air to actual mass (m) of 
body 

ae = reference area and reference length, respec- 
tively 

(1/2) pV? = dynamic pressure 

o = mL?2/I = nondimensional reciprocal of square of 
radius of gyration 

Cz, Cp = standard NACA coefficients for lift and 
drag forces perpendicular and parallel to 
flight path, nondimensionalized by divid- 
ing by [(1/2)p V2] 

Cis = standard NACA moment coefficient about 
center of gravity, nondimensionalized by 
dividing by [(1/2)pV?SL] 

Cy, Ce = standard NACA coefficients for forces nor- 


mal and parallel to principal body axis 
(x), nondimensionalized by dividing by 
[(1/2)p V2S] 
= derivative with respect to time (t) 
derivative with respect to nondimensional 
distance ~ traveled along flight path, see 
Eq. (14) 
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For Eqs. (22) to (24), when a = w/u — 0 and q = d6/dt > 0, 


X, = —puS(Cc), — Cc, > 0 

X, = —p(u/2)S(0C./da) = —p(u/2)S[(OCp/da) — C,] = 
— p(u2S/2) 

X, = = —p(uSL 2)Co, 

Z, = —puS(Cy), — p(u?/2) S(OCy/ou); Cy, ~ O 

Zy = —p(u/2)S(0Cy/0a) = —p(u/2)S[(OC,/da) + Cp] = 
— p(u?2S/2) 

Z, = = —p(u/2)SL Cz, 

Z; = —p(u?/2)S(0Cy/Ow) ~ O 

M, = puSL(C,,), + p(u?/2) SL(OC,,/Ou); = 0 

M,, = p(u?/2)SL (0C,,/dw) = (puSL/2)Cne 

M, = = p(u/2)SL? Ca 

M; = p(u?/2)SL(OC,,/Ow) ~ (pSL?/2) = 


(pSL?/2) Cme 


INTRODUCTION 


— MAIN PURPOSE of this paper is to discuss the 

range of validity of the usual linearized equations 
describing the dynamic longitudinal stability of a 
nonrolling ballistic missile that is attempting to follow 


Axes always tangential to the instantaneous flight path 
or actual trajectory. 


Fic.2. Axes fixed in body and coinciding with the principal axes 
of inertia. The y axis is orthogonal to and projecting out from 
the diagram. 
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a nearly zero angle of attack flight path while it is either 
accelerating or decelerating. 


ANALYSIS BASED ON SMALL PERTURBATIONS ABOUT 
FLIGHT PATH TRAJECTORY 


Both Friedrich and Dore! and Allen? developed their 
dynamic longitudinal stability equation by considering 
the zero thrust flight path trajectory equations which 
have the lift (C,) forces always normal to the flight 
path and the drag (Cp) forces tangent to it, as defined 
by (see Fig. 1 and Symbols) 


mV = —(1/2)pV2SCp — W sin y (1) 
mVy = (1/2)pV?SC_ — W cos y (2) 
+ & = (1/2)pV2SLC,, = 16 (3) 


The most fundamental assumption. then, is that 
W/pV2S = 0 (4) 


so that we may write, at least for hypersonic speeds and 
the nearly zero lift trajectory, 


(2mV /pV2S) ~ —Cp constant (5) 
(Qmy/pVS) ~ ~ (6) 


This means that we can only consider hypersonic speeds 
so that Cp and Cz, are approximately independent of 
flight speed or Mach Number. This same assumption 
also applies to C,,, but experience has shown that we 
must express it as®* 4 


Cn + Cnal(L/V)a] + (y +4) (7) 


As a direct consequence of these assumptions Eqs. (2) 
and (3) may then be combined to eliminate y and obtain 


a+ Da + Ma(t) = 0 (8S) 
= w(V/L) [Cra — + Cg) ] (9) 


M(t) = V/L)? BOC ng Cra | + 
\d[uCra(V/L)]/dt} (10) 


uw = (pSL/2m), o = (mL?/I) (11) 


We see that Eq. (8) is identical to Eq. (2) of reference 2, 
and also to Eq. (13) of reference 1 (if the additional 
term C,,, is added). It is important to note that terms 
written as C;, @ in reference 1 must be replaced by 
C,,, which cannot be a function of @ itself, otherwise 
the desired equation would not be linear in a@. 

It is obvious that Eq. (8) is not in the proper form for 
discussing the destabilizing effect of high rates of de- 
celeration, since the coefficients D and M depend di- 
rectly upon the rapidly varying velocity. 

In reference | the transformation 


a(t) = a(t) exp {| - f D(r) 2|art (12) 


is introduced, thereby reducing Eq. (8) to 
& + [M — (D/2)? — D/2]a(t) = 0 (13) 


which is now simpler to check for stability but again is 
not satisfactory for developing a general stability cri- 
terion, since the coefficient still depends directly upon 
the rapidly varying velocity. 

A better transformation for our purpose is the usual 
nondimensionalization that replaces the time by the non- 
dimensional variable £, which represents the number of 
body lengths traveled along the unperturbed trajectory 
or average flight path. This is more general than the 
transformation used in reference 2, since it is applicable 
to any arbitrarily curved trajectory 


t ) 
g = (1 1) f V(r)dr, (dé/dt) = [V(t)/L] | 
(*) = [d( )/dt] = (V/L) [d( )/dé] = rp (14) 
(V/L) | 
(V/L)? ( 4 (VV’ ‘L?) ( 
so that Eq. (8) is transformed into 
a” + ba’ + ca(t) = 0 (15) 


= — o(Cn, + + (V/V) (16) 


c(é) u[—oCn, Cr, + (V’ + 
(17) 


It then follows from Eq. (5) that we may write 
(V’/V) = L(V/V?) = —pCp (18) 


Now Eq. (15) is in the best form for developing a 
simple stability criterion, since 6 and c are independent 
of V(€) and depend only upon the monotonic variation 
of p(é) and the aerodynamic coefficients themselves. 
Eq. (15) is equivalent to Eq. (6) of reference 2, the 
main difference being that £ is the nondimensional dis- 
tance along the average curved flight path rather than 
the actual decrease in altitude. Allen® obtains approxi- 
mate solutions for his equation for a straight flight path 
by assuming an exponential variation of air density with 
altitude, which in our case would be written as 


= {1 — exp[—BL(E + &) sin y]} 
po ~ 0.0034 Ib. sec.? ft.-4, B = > (19) 
(1/22,000) 


The approximate solutions presented by Allen* show 
the direct effect upon the amplitude of subsequent 
oscillations of a because of a small initial error in @ 1f 
a is assumed to be zero. The effect of an initial error 
in both @ and @ can be evaluated directly from Eq. 
(15) by utilizing the inequality derived by Laitone.’ If 


= — — 6'(8)/2} = fnin > 0 (20) 


then it is shown ® that as long as f varies monotonically 
(that is, it is either never decreasing or never increasing, 
even though varying), the magnitude of a is bounded by 


= € 
la(é)| < V (f(0)a(0)? + }a’(0) + a(0) [b(0) fin exp [—b(é) /2|dé (21) 
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If the results of reference 5 were applied directly to 
Eq. (8) it would not be possible to have f(t) vary 
monotonically if the velocity were varying rapidly. 
However, {(é) can generally be shown to be monotonic 
when applied to Eq. (15). For example, if the aero- 
dynamic coefficients are assumed constant and p(£) is 
defined by Eq. (19), then the f(é) defined by Eqs. (15) 
and (20) is monotonic for a straight flight path. 


ANALYSIS BASED ON STABILITY AXES FIXED IN THE 
Movinc Bopy 


In order to evaluate the restrictions imposed by the 
approximations introduced in Eqs. (4) to (7), the pre- 
ceding derivation will now be repeated using Euler's 
classical dynamic equations of motion, based on axes 
fixed in the body so as to be always coincident with 
the principal axes of inertia. Consequently, the linear- 
ized forces now must be expressed as normal (Cy) and 
tangent (Cc) to the principal body axis (x), see Fig. 2, 
Symbols, and reference 3 or 4, 


m(u + wq) =~ + X,w + Xq 


— Weos w f gdt — (1/2)pV?SCe, (22) 
0 


m(w — uq) ~ Z,6u + Z,w + Z,q 
— W sin vw f gdt (23) 
0 
Ig ~ + M,w + + (24) 


for the nonrolling missile having only longitudinal 
oscillations. 

Since the forces and moments have been linearized 
by assuming that w and q are both small perturbations 
on the nearly zero lift flight path, we can linearize Eq. 
(22) simply by neglecting (wg). Then Eq. (22) is com- 
pletely independent of Eqs. (23) and (24) as long as Ce 
does not vary appreciably with small changes in u, w, or 
q; therefore, 


(25) 


using the same approximation as in Eq. (4). Similarly, 
corresponding to Eqs. (6) and (7), we have, to the same 
order of approximation, 


(2mi/pV?S) ~ —Cc, = constant 


2m(ug — = Cy,w (26) 
(21/pV*SL)q W + Cris Cmq[(L/u)q] (27) 


These linearized equations assume, as before, that the 
aerodynamic coefficients do not vary directly with the 
forward speed or Mach Number, so the analysis is again 
limited to hypersonic speeds (wu or V large). This is 
also necessary because of the introduction of Eq. (4). 

Then we can combine Eqs. (26) and (27) to eliminate 
q, and introduce Eq. (14) to eliminate the time, and 
finally obtain 


(28) 
(29) 


w” + bw’ + cw(t) = 0 


bi(€) 
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= wu(—oCy, — woCn Cy,,) + 


[d(uuCy,,) dé| (30) 
Now for small @ variation we can write 
w=au, w' =a'u+an’| 
” ” (31) 
we = au + 2a'u’ + an" | 
uCy, = UCn,, = Cn, (32) 


However, as will be discussed later, it entails an addi- 
tional approximation to write 


= (U?/L)Cns ~ Cn, = OCn/WaL/u) (33) 
since « is large and varying rapidly. Introducing 
Eqs. (31) to (33) into Eq. (28) results in 

a” + boa’ + ca(t) = 0 (34) 


= ulCv, — o(Cn, + + 2(u’/u) (35) 


c(t) = {ul—oCn, — +[d(uGy,)/dé]} + 
{ulCv, — o(Cn, + (u’/u) + (u"/u)} (36) 


It is seen that for small a Eq. (35) is identical to 
Eq. (16), since 


Cy, = Cr, + Co (37) 

while from Eq. (5) or Eq. (25) we have 
V/V = = L(a/u?) = — pCp (38) 
because Cp = Cc cosa + Cy sin al (39) 


= Co+ Cy, 0? = Co f 


However, Eq. (36) is not the same as Eq. (17) even 
for the simple case of rapid deceleration as defined by 
Eqs. (5) or (25). The difference being that Eq. (36) 
contains an additional term u*CpoC,,, since 


(u"/u) = (L/u)* — = 


— 


which may not be negligible for high rates of deceler- 
ation. 

If these higher order terms containing yu? and p’/p are 
neglected, then Eqs. (15) and (34) both reduce to 


a” + Cp + lee + 


u[—oCn, =0 (41) 


Now application of Eqs. (19) to (21) shows that as long 
as 


i> [—oCn,] > om Cp + Cn) ] > 0 
(42) 


then @ is bounded, and a damped stable oscillation will 
occur during rapid descent along a straight line tra- 
jectory through the atmosphere defined by Eq. (19). 

In these expressions the aerodynamic coefficients were 
considered constant, as would be the case at hypersonic 
speeds. However, Eqs. (17) and (36) both indicate 
that the first effect of a varying coefficient will appear 
as the variation of the lift curve slope with the nondi- 
mensional distance traveled, C he * However C,, could 
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DYNAMIC LONGITUDINAL 


depend upon & only if it varied with the velocity or 
Mach Number, 1. Consequently, since the speed is 
rapidly varying, Eqs. (6) and (26) must be modified in 
the following manner 


Cr = (0C_/Oa)ba + (0C,/0M) 6(V/a) = 
Ci, a + a) 


This would prohibit the reduction of the equations into 
an equation containing only a, even if the speed of 
sound (a) and the average speed are considered constant. 

Therefore, in Eq. (15) we have only justified the 
terms that appear in Eq. (41); that is, we must con- 
sider the aerodynamic coefficients constant and neglect 
the terms containing p* or p’. Fortunately, however, 
Allen® has already shown that numerically these very 
same terms are negligible in comparison to those re- 
tained in Eq. (41) for the ballistic missile at hypersonic 
speeds. 

However, for transonic speeds the variation of the 
aerodynamic coefficients with Mach Number must be 
included; consequently, an equation cannot be devel- 
oped for a alone. In addition, Eq. (4) cannot be used 
for transonic flow, so the gravitational component must 
also be considered. 

If it is justified to use Eq. (4) in hypersonic flow, this 
implies that for a rotationally symmetric nonrolling 
missile the yawing oscillations are also predicted by 
Eq. (41). However, in transonic flow the equations 
for the yawing stability will differ from those for the 
pitching (longitudinal) stability because of the difference 
in the gravitational components. 


COMPARISON FOR THE CASE OF CONSTANT 
ACCELERATION 


All of the previous analysis was primarily concerned 
with the case where a /pu? was approximately constant, 
as defined by Eqs. (5) or (25). However, the deriva- 
tions of Eqs. (15) or (34) are not altered in the linear- 
ized case, if one considers the special condition of con- 
stant acceleration (that is, % constant, so that 7 is 
always zero). The only change required is in Eq. (5) 
or Eq. (25), which now must have the necessary thrust 
added to maintain constant acceleration while the 
velocity increases, as in reference 6. As shown by 
Oswald,* these additional terms do not directly affect 
the homogeneous form of the a equation so that the 
left-hand side of Eq. (17) in reference 6 is identical to 
Eq. (34) if one considers the aerodynamic coefficients, 
the density, and the acceleration as being constant. 
Now, however, there is a serious discrepancy between 
Eqs. (15) and (34). This is perhaps better illustrated 
when the time is the independent variable, in which 
case Eqs. (14) and (34) yield 


m(Vy + 6V) [d(yo + by) /dt] = {( 2) p(Vo + + — W cos (yo + = 


mV [d(6y) /dt] = {(1/2)pVe2SCz, [ao + da + 2a0(6V/Vo)] — Wcos yo + (W sin y)iy} = 
(1 2) pVo?S[2C_,(6V/ Vo) ba 
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a bs& = 0 (43) 
b3(t) = p(u/L) [Cv — + + (44) 


c3(t) = L)? — Cy, 
(u/L) [d(uCy ,) dt) + u(u L) X 
[Cv, — + Cn,)] + (ti/u)} (45) 


This is again in complete agreement with the left-hand 
side of Eq. (16) of reference 6, but in serious disagree- 
ment with Eq. (8) for the special case of constant 
aerodynamic coefficients, constant density, and con- 
stant acceleration. Eq. (8) does not directly contain 
the acceleration in the damping term or the square of the 
acceleration in the equivalent spring constant. The 
negative sign attached to the latter term could lead to 
divergence for sufficiently large acceleration or decel- 
eration, a condition which is not predicted by either 
Eqs. (8) or (15). It would be necessary to apply the 
stability criterion given by Eqs. (20) and (21) in this 
case. 


COMPARISON FOR THE CASE OF CONSTANT VELOCITY 


Even for this simplest condition, consisting of small 
perturbations about a constant velocity flight path, it 
will be shown that there is a fundamental difference 
between the two methods of deriving the stability 
equations. There is an error involved that can lead 
to serious discrepancies with large drag bodies. For 
example, when the acceleration is negligible, Eqs. (15) 
and (34) become identical, except for the fact that Eq. 
(34) contains Cy, in place of the C,, occurring in Eq. 
(15). In view of Eq. (37), it is obvious that for a 
large drag body Eq. (34) will predict more damping and 
a larger spring constant than will Eq. (15). 


DISCUSSION AND CONCLUSIONS 


An analysis of the linearized solution for the case of 
constant average velocity along the flight path will now 
provide the important clues leading to an explanation 
of the sources of disagreement between the two meth- 
ods. The two fundamental causes of the different 
answers will be found to be the use of Eq. (4), that is, 
the neglect of the gravitational force component, and 
the variables selected in the linearization itself. 

When all the terms defined in the Symbols are substi- 
tuted into Eqs. (22) to (24), one obtains a stability 
quartic for the case of constant average flight speed 
and constant aerodynamic coefficients, as shown in 
references 3 or 4. The same determinant can be ob- 
tained from Eqs. (1) to (3) only if one uses Lagrange’s 
method of small perturbations, as shown in references 
7 or 8. For example, in our simple application where 
Crq is assumed to be independent of Mach Number, we 
still must write Eq. (2) as 


] (46) 
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Since Vo, a, yo are the constant average values along 
the flight path, then Eqs. (1) and (2) may be written 
for the average conditions as 


W sin yo = —(1/2)pVo?SCo, (47) 
W cos yo = (1/2)pVo?SCi, = (48) 


Consequently, even though C, is not considered a direct 
function of the variation in forward speed, still 6V’/ Vy 
appears in each of the linearized equations and thereby 
prevents their reduction to an expression containing 
a alone, unless C;, = 0, in which case the term Cp,67 
still prohibits the solution for @ alone unless, in addi- 
tion, Cp, = 0. Obviously, this cannot be assumed 
for a blunt body. 

Now, however, Eq. (46) reveals the important fact 
that Cp) will not combine with C,, as a term in the 
stability quartic if the variables are taken as (I, a, y). 
The combination found in Eqs. (35) and (37) will occur 
only if the variables are changed to (V, a, @), where 9 is 
defined by 


6=yta= 
6 = yta=q 
that is, only if the pitching velocity replaces the vari- 
ation in flight path direction as a new variable. This 
is immediately evident when Eq. (49) is combined with 
Eq. (46) to yield 
[d(6y) dt| = [d(60) dt| — [d(éa) = 
(pVoS/2m) [2C,,(6V) Vo) + 
(Cra + Co,)ba — Cp,60] 


Therefore, one derives the same linearized stability 
quartic only if one uses (V, a, #), as in references 7 or 8, 
or (uv, w, g) asin references 3 or4. However, the linear- 
ized stability quartic derived from (I, a, y), although 
differing in the manner discussed from the classical one 
given in all the references, is still valid, and can be just 
as applicable in some cases. For example, when 
Cra > Cop,, then, obviously, dy, the instantaneous flight 
path, and the actual motion of the center of gravity 
could deviate considerably from the average flight path, 
especially for the long period or phugoid motion, where 
da — 0. Although in this case 6y ~ 66, still Eqs. (15) 
and (34) would be equivalent, since Cy, ~ Cry > Cop,. 
Similarly, in the limiting case, when C,;, — 0, then 
5y — 0, so that both terms containing Cp, from the 
gravitational component should be neglected sirce da ~ 
60. Physically, this corresponds to the limiting case 
of the nonlifting blunt body which is oscillating about 
its unperturbed center of gravity, as would be the situ- 
ation in the wind-tunnel test of any freely oscillating 
model. It is obvious that Cp, could not contribute to 
the damping if the center of gravity did not oscillate nor- 
mal to the average wind direction. Consequently, 


(49) 


(50) 


whenever 67 = 0 = C,,, Eq (15) is more suitable than 
Eq. (34), especially for a blunt body having a large drag. 

It is now evident that the necessity of neglecting the 
gravitational component in order to obtain a linear 
expression in @ alone has resulted in the unavoidable 
Eq. (15) has neglected 


omission of various terms. 
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Cp,6y so it would be more applicable whenever 6y ~ 0 
and 6a ~ 60; that is, whenever C,, — 0, or the aero- 
dynamic moments are much greater than the lift forces. 
On the other hand, Eq. (34) has neglected Cp,60; conse- 
quently, it is best suited whenever 60 ~ 0 and é6y =~ 
—éa; that is, whenever the aerodynamic moments are 
small in comparison either to the moment of inertia, 
or to the lift and drag forces. 

These conclusions were obtained for nearly constant 
velocity. For large rates of acceleration additional 
difficulties are introduced because w and & are no longer 
directly related as shown by Eq. (31). For example, a 
comparison of Eqs. (28) and (34) shows an entirely 
different effect of acceleration upon w and a, the two 
equations become identical only when wu is constant. 
A somewhat analogous situation is discussed by Reed*® 
who shows that the damping of variables such as a, 4, 
or & is different, and each is differently affected when 
the oscillation frequency varies for any system that 
has its coefficients varying with time. 

Fortunately, one of the effects of acceleration— 
namely, the difference between C,,, and C,,.—-vanishes 
at least for a blunt heavy body since Cnn, is approxi- 
mately zero in hypersonic flow. However, for a very 
light body, Cy, and C,,. are directly related to the 
apparent mass of the surrounding air; therefore, in this 
case, Ww, and not a, is the fundamental variable which 
is best suited for high rates of acceleration. On the 
other hand, both methods justify the application of Eq. 
(41) to a large drag blunt body having a high rate of 
deceleration at hypersonic speeds since it is the first- 
order solution by either method, and the order of mag- 
nitude of the terms retained has been justified by Allen.? 
In this case Eq. (42) provides the appropriate stability 
criterion. 
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Supersonic Airfoil Performance With Small 


Heat Addition’ 


ARTUR MAGER* 
Marquardt Aircraft Company 


SUMMARY 


An analytical method is presented which permits a very rapid 
evaluation of the aerodynamic effects arising from the addition 
of small amounts of heat near supersonic two-dimensional air- 
foils. This method applies to shockless inviscid flow without heat 
conduction. Also, the mechanism by which the desired heat 
addition is achieved is not considered. 

It is shown that even small amounts of heat generate a sub- 
stantial pressure rise and thus cause appreciable changes in the 
aerodynamic coefficients. The results of this analysis compare 
favorably with those obtained by a more accurate, but also more 
tedious, graphical method of characteristics. 

Two possible modes of application to an airplane design are 
considered from the energy requirements standpoint. In this 
connection, it is shown that the decrease of the required wing 
area resulting from heat addition may, in some cases, lead to 
savings in the rate of the fuel consumption. In general, however, 
one should not expect any substantial reduction in energy re- 
quirements resulting from the application of the wing heat ad- 


dition. 
SYMBOLS 

a = speed of sound 

Cc = coefficient 

= chord 

D = drag 

E = engine performance function 

G = heat addition function (sec./ft.)? 

g = gravitational acceleration 

H = total rate of heat addition per unit of span, 
B.t.u./(ft.sec.) 

h = projection on y axis of that part of forward running 
Mach wave which is inside the heat addition re- 
gion 

h; = heat addition region height 

I = unit step function 

J = ft.-lb. equivalent of B.t.u. 

L = lift 

M = Mach Number 

p = static pressure 

Q = strength of a unit-long heat source line, 
(ft.lb.)/(ft.*sec. ) 

R = gas constant 

s.f.c. = engine specific fuel consumption 

= total temperature 

t = static temperature 

u,v = velocities in x and y directions 

W = airplane weight 

w = mass flow per unit span 


Presented at the Propulsion Aerodynamics Session, 26th An- 
nual Meeting, IAS, New York, Jan. 27-30, 1958. Revised and 
received September 5, 1958. 

+ This investigation has been sponsored by the Marquardt 
Aircraft Co. 

* Technical Consultant, Astro Division. 

The author wishes to express his thanks and appreciation to 
R, E. Fisher and Dr. C. A. Lindley for the many constructive 
suggestions offered him in the course of this work. 
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X = dimensionless distance along chord (x, 
x,y, 2 = coordinates 
a = angle of attack 


= Mach wave angle 
Y = ratio of specific heats 
A = increment due to heat addition 
6 = Dirac delta function 
€ = relative efficiency of wing heat addition 
¢ = angle between generated force and vertical 
n = ratio of fuel consumptions with and without wing 
heat addition 


rs) = airfoil and heat addition geometry factor 

K = ratio of wing frictional drag to total wing drag 

r = ratio of body to airplane drag 

p = density 

=c = 

o = angle between airfoil surface and chord 

T = airfoil thickness, relative to chord 

re) = fuel heating value 

y = ratio of body plus frictional drag to airplane drag 


Subscripts and Superscripts 


A = airplane 

B = body 

b = basic (at a = 0) 

D = drag 

f = at end of heat addition region 

H = heat addition 

1 = at start of heat addition region 

iL = lift 

0 = unperturbed, ahead of the airfoil 
= surface 

W = wing 

a = with varying a 

T = friction 
= at end of region I 

2 = at end of region II 

vd = break-even value 


bar above value signifies heat addition 


INTRODUCTION 


Fg the basic equations describing super- 
sonic flow with heat addition, also known as 
diabatic flow, have been established and studied for 
some time (primarily by B. Hicks and his coworkers) ,!~* 
only a few problems of practical interest have so far 
been considered. One reason for the relative inactivity 
in this field is the fact that the heat addition, with the 
resulting variation of total temperature and total pres- 
sure, complicates considerably the task of finding the 
appropriate flow field around even the simplest aero- 
dynamic configurations. Some progress towards sim- 
plifying this problem has been made by Pinkel and Sera- 
fini,» who have devised a graphical method of charac- 
teristics applicable to heat additions with continuous 
variation of total temperature in the streamwise direc- 
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Fic. 1. Supersonic line heat source. 


tion. Using this method, Pinkel and his coworkers 
have evaluated the pressure distribution and the ap- 
propriate aerodynamic coefficients for a simple bi- 
convex supersonic profile.6 However, in spite of the 
simplifications and graphical aids to computation in- 
cluded in reference 5, the amount of work required to 
carry out such calculations remained quite appreciable. 
Furthermore, once such work is performed for a given 
aerodynamic configuration, it is difficult to draw from 
it any general conclusions applicable to other aero- 
dynamic shapes, because a clear analytical relation be- 
tween the aerodynamic and the heat addition param- 
eters is impossible to establish. 

On the other hand, it may be expected that only 
small amounts of heat will ever be added externally to 
the stream (such was also the case in all but one of the 
examples of reference 6) so that a small perturbation 
procedure based on this fact should be applicable. In 
principle, such a procedure, while not as accurate as 
other more lengthy methods, leads to a certain essen- 
tial sorting process between the cause and effect, so 
that it should permit a considerably more direct inter- 
pretation and description of the external heat addition 
problem. In this connection, Tsien and _ Beilock’ 
presented results describing the effect of a linearized 
heat source on the uniform flow, and these, when prop- 
erly used, should be directly adaptable to such an ap- 
proach to the problem. 

In line with the above comments, the objectives of 
the present investigation are twofold. First, the heat 
source line concept of reference 7 is used in a linearized 
method for treating the heat additicn to shockless 
supersonic flow near an airfoil profile. This permits 
the establishment of the desired analytical relations 
between the changes in the aerodynamic coefficients 
and the amount and mode of the heat addition. Sec- 
ond, the relations so developed are used to indicate the 
overall efficiency of the process by evaluating the air- 
plane energy requirements ratio with and without the 
external heat addition. Quite intentionally, the pres- 
ent investigation concentrates on the examples and 
follows the general lines of the analysis of reference 6. 
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This permits not only the evaluation of the accuracy of 
the linearized method, but also some extension, clari- 
fication, and generalization of the results. 

It should be stressed that the method here proposed 
should only be used for a preliminary, but conveniently 
rapid, appraisal of the problem. It applies to shockless 
inviscid flow, with small disturbances, in which the 
heat conduction does not occur. Also, the actual 
means by which the desired heat addition will be ac- 
complished and possible restrictions arising therefrom 
are not considered.* 


LINEARIZED EFFECT OF HEAT ADDITION 
NEAR SUPERSONIC AIRFOILS 


Pressure Distribution Change 


As mentioned in the Introduction, Tsien and Beilock’ 
give the expressions for the effect of a line source lo- 
cated in the origin of a uniform flow field, under the 
assumption that the total amount of heat added is 
small and the equations can thus be linearized. The 
appropriate perturbation quantities may be written as 


Gay? 
v= Gao? V 
p= (1) 


GyMobo{ 1 — Me? — X 
[5(y) — yVW Mo? — 


with G defined by 


> 
II 


G = I(y — 1)O/2ya0*poV Me? — 1] X 
+ yV Me? 1) 


where 6 is the Dirac delta function, and /() is the unit 
step function. From expressions (1) one can see that 
the disturbances due to the heat source are thus limited 
to a very narrow region near the lines y = x iV Me — 1; 
that is, the Mach lines of the basic flow. The density, 
of course, is also affected by the heat wake extending 
along the x axis. All this is illustrated schematically 
in Fig. 1. Furthermore, examination of reference 7 


* In the interval between the completion of this work and its 
presentation, the author’s attention was called to two additional 
papers dealing with the same subject. One of these, by C. 
Gazley, Jr., ‘Linearized Solution for Heat Addition at the Sur- 
face of a Supersonic, Airfoil’? (RAND Corp. RM-1892, ASTIA 
AD 133025, November, 1956), employs a “‘piston’’ concept to 
translate the effect of heat addition into an effective deflection of 
the flow. The second, by W. W. Willmarth, ‘‘The Production 
of Aerodynamic Forces by Heat Addition on External Surfaces of 
Aircraft” (RAND Corp. RM-2078, ASTIA AD 150681, De- 
cember, 1957), makes use of the analogy between the heat source 
and the fluid-mass source. Since both of these reports resort to 
linearized equations, they arrive at essentially the same results 
as those presented here. 

In addition, Willmarth also formulates a simple model of heat 
addition large enough to result ina shock. The pressure rise per 
amount of heat added is, of course, in such a case lower than that 
obtainable from the linearized analysis. Therefore, from the 
overall efficiency standpoint, one may expect the linearized 
analysis to give an upper bound on the values which will be 
actually realized. 


dise 


easi 
wit 
met 
sour 
add 
turb 
to 
the 
zero 
one 
of a 
seco 
cally 
bout 
ity 
iden 
side 
heat 
ever 
tion 
com 
poin 
Eqs. 
tive 
W 
will | 
from 
lines 
warc 
num 
Mac 
the 1 
axis 
regic 
Obse 
inter 
to t 
This 


must 


100 
| v is 
i 
| 
Cas 
B 
D 
a 
E 
F 


SUPERSONIC AIRFOIL PERFORMANCE 101 


discloses that wu, p, and p are even functions of y, while 
7 is odd—that is, 


v(x, vy) = —v(x, — y) (2) 


These preliminary, but basic, results can now be 
easily applied to the problem of supersonic airfoils 
with heat addition. Naturally, since this is a linearized 
methcd—that is, since the effects of individual heat 
source lines, angle of attack, thickness, camber, etc., are 
added arithmetically—to compute the pressure per- 
turbations due to the heat addition one does not need 
to consider an actual airfoil; it will suffice to visualize 
the heat addition region as distributed over a chord at 
zero angle of attack. To be consistent with reference 6, 
one thus imagines a rectangular heat addition region 
of area hi(xy — x;,) delivering B.t.u.’s of heat per 
second and per unit of span length. This is schemati- 
cally illustrated in Fig. 2. Furthermore, since the 
boundary conditions require that the perturbation veloc- 
ity v vanish on the chord, one must also imagine an 
identical heat addition region, located on the opposite 
side of the chord, and delivering the same amount of 
heat. For practical purposes, of course, since p is an 
even function of y, one only has to double the perturba- 
tion obtained from the actual heat addition region. The 
computation of this pressure perturbation acting at any 
point on the chord may be carried out directly from 
Eqs. (1), after the Q used there is replaced by an effec- 
tive value Q.¢¢,, which will now be determined. 

We take c as a unit of length, so that each source line 
will have the strength JHc/[(x,;— x,)h;]. It is obvious 
from the previous remarks that the only heat source 
lines affecting a given point are those lying on the for- 
ward running Mach waves from that point. Their 
number will clearly be given by the length of such 
Mach waves inside the heat addition region, divided by 
the unit length c. Thus, if / is the projection on the y 
axis of the Mach wavelength inside the heat addition 
region, then this number may be written as //(c sin 8). 
Observe now that the elemental area cdx, which is of 
interest here, is larger by 1/sin 8 than the area normal 
to the Mach waves, which is implied in Eqs. (1). 
This means that the effective pressure perturbation 
must be correspondingly decreased by sin 6. Conse- 


quently, taking account of all these effects, we write 
Qett. = — (h/c sin B) sin B (3) 


Upon the substitution of this value in Eqs. (1), the pres- 
sure perturbation at any point on the chord (or on the 
airioil surface, within the assumptions of the linearized 
theory) is given nondimensionally as 


Ps/Po = {(y 1) MocJ Cy /|(xy — Xz) x 
V ygR(Mi? — 1)]} (4) 


Here Cy, the coefficient of heat addition from reference 
6, is defined by 


H = — 1] 


Incidentally, it should be obvious irom the above der- 
ivation that Eq. (4) is also valid for any other heat 
addition region of arbitrary shape and/or Jocation, as 
long as one substitutes the new cross-sectional area of 
such a region for the rectangular area (x, — x,)hy. 

One thus sees that the pressure change due to the 
rectangularly shaped heat addition is intrinsically 
dependent on the ratio (h/h,;) and, accordingly, recog- 
nizes three distinct regions in which this ratio has a 
certain form. This is illustrated schematically in Fig. 
3. In region I, extending from x; to x, with x, = 
x, + hy Me 1, the surface perturbation pressure 
increases linearly with the distance because h/h; = 
(x — x;)/(%1 — x). In region II, extending from x, 
to x,y, the surface perturbation pressure is greatest, and 
constant, since h = h;. Finally, in region III, extend- 
ing from x, to X2, with v2 = x, + hiv Me — |, the sur- 
face perturbation pressure decreases linearly as h/h; = 
1 — (x — xy)/(x2 — Xy). 

With these variations of (/h,), one can easily com- 
pute the surface pressure perturbation from Eq. (4). 
This has been done for the examples of reference 6 
concerning a 5 per cent thick airfoil, with all the per- 
tinent data listed in Table 1. The representative re- 
sults of such computations are compared to those ob- 
tained by the more accurate method of characteristics® 
in Figs. 4(a) and 4(b). As may be seen from Fig. 4(a), 
the agreement is quite good when the assumptions of 
the linearized theory hold—that is, when the heat 


TABLE 1 
T;/To H/H Q/Q 
Case Mo xi/C x7 /C h;/¢ (a = 2°) Cu (Case B) (Case B) Remarks 
A 3 0 0 0 1.0 0 0 0 No heat addition 
B 3 0.356 0.690 0.07 1.126 0.0201 1.0 1.0 All generated pressure 
felt by the airfoil 
 & 3 0.356 1.0 0.07 1.243 0.0887 1.925 1.0 Some generated pressure 
lost 
D 5 0 0 0 1.0 0 0 0 No heat addition 
E 5 0.356 1.0 0.07 1.243 0.166 8.260 4.282 Some generated pressure 
lost 
F 3 0.356 1.0 0.085 1.243 0.0190 0.945 0.982 Some generated pressure 
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addition coefficient is small. For large values of Cy, 
Fig. 4(b) indicates that the linearized theory tends to 
overestimate the pressure rise due to the heat addition. 


Lift and Drag Coefficient Change 


The results of the previous section can now be used 
to compute the changes in the lift and drag coefficients. 
For the sake of clarity it is convenient to assume here 
that the airfoil is at zero angle of attack and to deter- 
mine this effect later. The pressure perturbation 
due to the heat addition generates a force increment 
acting locally on an elemental surface inclined at an 
angle o with respect to the chord. It is trivial to show 
that the elemental increase of the lift coefficient is 
independent of the airfoil shape and can thus be written 
as 


= [(2ps)/(yMe%po) dX (6) 


Eq. (6) must now be integrated, taking into account 
the proper variation of (/h;) in each region. When 
this is done, the total change of the lift coefficient due 
to the heat addition becomes 


It should be noted that Eq. (7) pertains to the case 
where all the pressure generated by the heat addition is 
felt by the airfoil—that is, when X, < 1.0. It may 
sometimes happen that this will not occur, and part, or 
all, of region III falls outside. By proper change of inte- 
gration limits, it is easy to show that, for XY, = 1.0, 
the total change of the lift coefficient becomes 


ACL. 
ACz. — [0.5(h,/)V Me — 1/0 — (8) 


and consequently, as might have been inferred from 
basic principles, such a location of the heat addition 
region results in a loss of the realizable lift increase. 

As is well known, the elemental increase of the drag 
coefficient is a function of the airfoil shape, and, for the 
linearized theory, it is given by 


dCp [(20ps) /(yMo"po) (9) 


In general, therefore, the integration of the total drag 
increase due to the heat addition cannot be performed 
unless the airfoil shape is first assumed. However, 
for the circular are biconvex airfoil, such as used in 
reference 6, one can perform the integrations without 
much difficulty. For such an airfoil, if the thickness is 
small (e.g., 7 < 1.0), it is elementary to show that 


o = —27(1 — 2X) (10) 
and, consequently, Eq. (9), with proper variation of 
(p,/po) and with X» < 1.0, integrates to 
= —27(X; + X2 — — OAC,» (11) 


with 0 being of O(7) and signifying a certain combined 
airfoil and heat addition region geometry factor. 
Thus, it is obvious that the greatest reduction in drag 
will be obtained when X; and X, are large—that is, 


when the heat addition region is near the trailing edge 
of the airfoil. 

Of course, when region III falls outside the airfoil, as 
is the case in some examples of reference 6, a change in 
the integration limits must be introduced and then the 
appropriate expression for the total drag decrease 
becomes 


ACp,»)xyar = —274{X; 
[(hi/W Me — 1/11 — + 
(hi/3c)V Me? — 1 — (12) 


Once the lift and drag coefficient changes at zero 
angle of attack are known, the effect of angle of attack 
can easily be evaluated. All one has to do is to imag- 
ine that the heat addition results in a force inclined at 
angle ¢ to the vertical. As the angle of attack is 
changed, since the heat addition region remains fixed 
with respect to the airfoil and since in the linearized 
theory the effects are separately additive, this force, 
unchanged in magnitude, simply rotates with the air- 
foil. If the angle of attack is taken positive clock- 
wise, then 

¢=% — a= —[(ACp/AC,), + a] (13) 


and the elementary resolution into components shows 
that 


AC, = (cos ¢{/cos {))ACz, » = ACz.» (14) 


ACp = (sin ¢/sin = ACp,» X 
[1 + a(AC; ‘AC p)o] (15) 


One thus sees that, since ACp,, < 0, the angle of 
attack has a pronounced influence on the decrease of 
drag due to the heat addition and that the greatest 
reduction of the drag coefficient for this airfoil will be 
obtained when a is negative. In general, there will be 
some angle of attack at which a reversal in the direction 
of the drag change will occur, and this angle, for the 
present example with X. < 1.0, is given from Eg. (11) 
and (15) as equal to 0, because 


ACp —AC,,(9 (16) 
Since the expressions for the lift and drag coefficients 
without the heat addition are well known and can be 
written as 
C, = 4a/VM? — 1 (17) 
Cp = A[a? + Me? — 1 (18) 
then the total values, with heat addition, can also be 
obtained. For X2 < 1.0, these are 
C, = (4/V Me — 
— V ygR]} (19) 
(4/V My? — 1){ a? + (4/3)7? — 
[((y — I) ygR(O — a) (20) 


Cp 


It should be noted that the variation in angle of 
attack causes changes in the total temperature which 
will finally be reached, even if the heat addition region 
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and the total amount of heat added remain unchanged. 
This occurs, of course, because of the variation of 
mass flow passing through the heat addition region. 
For small height of the heat addition region h;, these 
changes can easily be estimated by assuming a one- 
dimensional] flow at proper surface inclination and lin- 
earizing the resultant equations. Under these assump- 
tions, the final total temperature reached can be ob- 
tained approximately from 


T,/To = 1 + [(T;/To)» — 1] X 
(i ~ Me 1) 


Using the above results, the cases listed in Table 1 
were evaluated for C, and Cp and then compared to 
their corresponding values given in reference 6. This 
comparison is shown in Figs. 5(a) and 5(b), and it may 
be seen from these Figures that the agreement between 
the linearized theory and the method of character- 
istics is quite good. 

In Fig. 6, the effect of the mode of the heat addition 
on the aerodynamic performance at M/) = 3 and 7 = 
0.05 is shown by plotting the results for cases B and C. 
It will be noted from Table 1 that case C, although 
having the same rate of heat addition as case B, differs 
from case B because it has nearly twice as much total 
heat added. Furthermore, while in case B all the gen- 
erated pressure is felt by the airfoil, in case C some of 
the generated pressure is lost. Corresponding to these 
conditions, it may be seen from Fig. 6—although case 
C shows slightly greater changes in lift and drag co- 
efficients than case B—that this increase is not pro- 
portional to the increase in total heat added. As 
mentioned previously, since the changes in lit coeffi- 
cient are unaffected by the angle of attack and since 
the reduction of the drag coefficient tends to decrease 
with angle of attack, one may expect that the maxi- 
mum lift-drag ratio with heat addition wil] tend to 
occur at lower angles of attack than that without the 
heat addition. This trend is shown in Fig. 6, where 
a for best Z/D has decreased for cases B and C by 
approximately one degree from its value for case A, 
when there is no heat addition. 

As may be seen from Figs. 5 and 6, none of the ex- 
amples used had heat addition high enough to result 
in actual thrust. It may be interesting to inquire 
what is the necessary value and mode of heat addition 
to give sucha condition. By considering for the sake of 
simplicity only the case when a = 0, one finds, from 
Eq. (20), that the necessary requirements for the bi- 
convex circular arc profile can be written as 


Cu(Xi+ X2— 1) > 
(4 V ygR)/[BJ(y — 1)] (22) 


Consequently, for higher flight Mach Numbers, more 
heat will be required to obtain net thrust from the ex- 
ternal heat addition to the supersonic airfoils. It 
should be realized, of course, that the requirement indi- 
cated by relation (22) must be applied cautiously be- 
cause of the linear character of the analysis and be- 
cause of the assumed shockless and inviscid flow. 
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COMPARISON OF THE ENERGY REQUIREMENTS 


General Assumptions 


In order to estimate the possible limits on the appli- 
cation of the external heat addition to an airplane as 
a whole, one must consider its cost in terms of the 
necessary energy requirements. As pointed out in 
reference 6, such considerations can be carried out in 
two different ways. One of these is to assume that the 
airplane with heat addition flies at the same angle of 
attack but requires a smaller wing area to carry the 
same load as the airplane without the heat addition. 
This decrease in wing area results in a corresponding 
decrease of drag and, thus, causes a change in the energy 
expended by the engines for flight. Of course, the 
energy used in the external heat addition process must 
also be taken into account. Another way is to assume 
that the wing area of the airplane with the external 
heat addition remains unchanged, but now, with heat 
addition, the same load can be carried at the appro- 
priately lower angle of attack. Again the energy re- 
quirements are changed, and these may be better, or 
worse, than for an airplane without the heat addition. 
Both of these ways will be considered separately. 

In carrying out the analysis, it was assumed, just as 
in reference 6, that the hypothetical airplane has rec- 
tangular constant profile wings, unaffected by tip 
effects, and that no weight penalty is involved in the 
application of the externa] heat addition. The present 
analysis differs from that in reference 6, however, not 
only because it is now possible to use the previously 
derived analytical relations, but also because it assumes 
that the frictional drag coefficient, and not the total 
drag coefficient of the wing, is unchanged by the heat 
addition. It will be seen that this last difference is 
quite crucial, because it results in a certain charac- 
teristic ratio of frictional to total wing drag coefficients, 
which, in fact, determines the break-even point of the 
energy requirements. Because of this essential differ- 
ence in the assumptions used, something should be said 
about their relative validity. Unfortunately, however, 
at present, in the absence of the experimental data and 
in the absence ef boundary-layer solutions with ex- 
ternal heat addition, one is not able to evaluate whether 
this difference in the assumptions is really justifiable. 
Still, it should be noted that the present set of assump- 
tions is somewhat more consistent with the linearized 
analysis. This is so because the heat addition, when 
properly applied, will cause a pressure rise as well as a 
pressure drop, thus somewhat canceling the possible 
changes in surface shear. Furthermore, the amount of 
heat added is small, thus minimizing the possibility of 
shock formation and avoiding the problem of shock- 
boundary-layer interaction. 

The energy requirements of an airplane are con- 
veniently given by the required fuel consumption rate. 
Consequently, by comparing this rate for an airplane 
without the wing heat addition, when all the fuel is 
used by the engines, and that with the external heat 
addition, when the fuel is divided between the engines 
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and the wing heat addition mechanism, one can obtain 
a measure of the relative energy requirements. 


Airplane Weight and Angle of Attack Invariant 

By slightly modifying the expression derived in 
reference 6, the ratio of the fuel consumptions with and 
without the external heat addition for the same weight 
airplane can be written as 


n = (L/D)a[(D/L)4 + (ECu/Mo?C,)] (23) 


where /£, the engine performance function, is nearly 
independent of Mach Number, altitude, and fuel used, 

E = 7200Vt,/yo(s.f.c.) (24) 
The total airplane drag, without the heat addition, 
consists of both the wing drag and the body drag. The 
body drag can be stated as a fraction of the total air- 
plane drag—e.g., Dg = \D4—while the wing drag can 
be divided into that due to friction and that due to 
pressure ; 

Dw = + Co, 


As a consequence of such a representation, one has 


= (1 r)] [(Cp + Co. = 
W=C/(1 — A)C, (25) 


With the external heat addition, at the same angle 
of attack, the body drag and the wing frictional drag 
coefficient are assumed to remain unchanged, while the 
pressure drag and the lift coefficients are affected. 
With these assumptions one has 


Dg and Dy = (W, + Cp, ,) 


In view of these modifications, the total airplane drag, 
with the external heat addition, becomes 


Ds = Wi [Co + Co, + 
[\/(1 — A)(EC/Cz)} (26) 
Since the required lift of the airplane must remain un- 
changed—that is, since Ly = L4 = W—Egq. (23) be- 
comes, after substitution of Eqs. (25) and (26), and 
some manipulation, 
n= [(1 — X 
[Co + + (ECu/Mo*)) (28a) 
Let « define the ratio of the frictional to total drag 
coefficients for the unheated profile ; 
k=Cp,,/=C_ or (equivalently) 

=Cp/(1 — «) (27) 
then, because Cp = Cp + ACp and C, = C, + ACz, 
one can write Eq. (23a), after substitution and rearrange- 
ment, as 
(n — »+)/(1 — A) = — + AC,)] X 

— «)] + ACp + (ECu/Mo*)}  (23b) 

The break-even point in fuel consumption rates 

occurs when 7 = 1. If the value of x, which gives this 

value at any given a, is designated as x*, then from Eq. 
(23b) one obtains 
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g* = — ‘(1 + €p) |} (28) 
where ¢€, and e, are defined by 
=  €p = (29) 


and are a measure of the relative efficiency with which 
the energy is converted into the aerodynamic forces by 
the external heat addition process and by the engines. 

One observes, from Eq. (28), that the break-even 
condition is completely independent of body drag. 
By solving Eq. (28) for (Cp/C,) and substituting into 
Eq. (23b), one can also express the ratio of the fuel 
consumption rates as 


n= 1— {[(« — «*)/(1 — «*)] X 
(1 — A) [1 + (C,/AC,)]}} (23ce) 


Consequently, since \ and « must be smaller than unity, 
the airplane with the external heat addition will show 
advantage in fuel consumption rate when «x is larger 
than x*. One thus desires to have «* as small as pos- 
sible. Keeping this in mind, one can now substitute 
the appropriate expressions for (Cp/Cz), ez, and ep, 
so as to get a better idea of the factors affecting x*. 
Assuming that ¥. < 1.0, this gives for the biconvex 
circular are airfoil, upon retention of second-order 
terms, 


= 2(y — 1)J/yEV ygR[1 — (1/Me2)] (29a) 


= 1 — (e,/a)[a? + (4/3)r?] X 
[1 + «(8 — a)] (28a) 


Since the terms involving 0 in Eq. (28a) are of second 
order, it is clear that x* is only slightly affected by the 
changes in pressure drag, and, also, if Mo? > 1.0, «* 
is practically unaffected by flight Mach Number. For 
a given e,, low values of «* will occur when angle of 
attack is small and when airfoil thickness is large. It 
should be noted that, at angles of attack for zero lift 
(a = 0), x*— — ©, thus indicating that any airfoil will 
show an advantage there. This, of course, results from 
the fact that any change in lift for such a gives an in- 
finite change of wing area and, thus, an infinite reduc- 
tion of energy requirements. 

Eqs. (28a) and (29a) are valid when the heat addition 
region is located properly on the airfoil and all the pres- 
sure generated by the wing heat addition is converted 
into lift and drag changes. In some of the examples of 
reference 6 this was not the case, and it may be inter- 
esting to see how x«* is affected. This is illustrated in 
Fig. 7 for cases B, C, and E. It may be seen from this 
Figure that case B, which satisfies the condition X2 < 
1.0, also has a greater range of a values in which the 
advantage in fuel consumption occurs. It can also 
be seen from this Figure that cases C and E, which 
are at two different Mach Numbers but have a similar 
mode of heat addition, show about the same range of 
advantage in fuel consumption rates. 

By substituting the appropriate values at each a into 
Eq. (23c), one can also evaluate the fuel consumption 
ratio 7. This is shown for cases B and C in Fig. 8. 
In computing this Figure, E corresponding to a typical 
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a DEGREES 
CASE E 
ADVANTAGE\_ 
DISADVANTAGE 
« DEGREES 
Fic. 7. Effect of wing frictional drag on advantage in fuel 
consumption. 


ramjet flying above 35,000 ft. and X = 0.25 were used. 
As can be seen from this Figure, not only does the break- 
even condition occur for a larger value of a for case B 
than that for case C, but also the disadvantage in fuel] 
consumption rate, when it occurs, is small, appreciably 
smaller than that for case C. 

Finally, the effect of airplane body drag on 7 is shown 
in Fig. 9. As can be seen from this Figure, the fuel 
consumption rate ratio is not very sensitive to the cor- 
rect value of body drag, particularly for angles of attack 
not too different from the break-even value. This fact 
emphasizes the importance of « in the evaluation of 
advantages to be gained from the external heat addition 
process. 

It should be pointed out that, for the above examples, 
the advantages in fuel consumption rate tend to dis- 
appear for angles of attack greater than about 3 to 3.5 
degrees. These values of a are close to the values for 
the best lift-drag ratio for an airfoil without the heat 
addition, and it may be expected that an actual air- 
plane will incorporate just these values for reasons of 
best economy. But in the above analysis the angle 
of attack is invariant, and, thus, the possible advantages 
in fuel consumption rate due to the external heat ad- 
dition becomes very marginal. 


Airplane Weight and Wing Area Invariant 


If the airplane wing area remains unchanged, then, to 
carry the same weight with wing heat addition, the lift 
coefficient must also remain unchanged and C, = c. 
The airplane can thus fly at a smaller angle of attack 
when the heat addition is used. 

Since the wing frictional drag was assumed to be in- 
dependent of heat addition, to be consistent one should 
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also assume that this drag component is independent 
of the angle of attack. Furthermore, in view of the 
invariant wing area, it is convenient to combine the 
wing frictional drag with body drag and express the sum 
as a fraction of airplane drag. Thus, one has 
+ Dz = Dz VD, 
By applying these changes, one can now write the air- 
plane drag without heat addition by modifying Eq. 
(25) as 
Ds = — ¥))(Co/Cr) (30) 


and with wing heat addition, from Eq. (26), as 


Ds = + — (31) 


But L4 = Ly = W, and the fuel consumption ratio is 
still given correct!y by Eq. (23). Consequently, upon 
substitution of Eqs. (30) and (31) into Eq. (23) one 
obtains, after some manipulation, 


n= 1+ — ¥) (1 — (32) 


where ¢, is again a measure of the relative efficiency 
with which the external heat addition process and the 
engines convert the supplied energy into the aero- 
dynamic forces, 


ia = [(Co — Cp) Mo?]/(ECu) (33) 


It is quite obvious from Eq. (32) that the break- 
even condition occurs when ¢, = ¢, = 1.0. Further- 
more, since Y must be smaller than unity, to obtain an 
advantage in the fuel consumption rate with wing heat 
addition, one desires values of ¢«, larger than unity. 
With these facts in mind, it is appropriate to evaluate 
the parameter €,. Taking again the biconvex circular 
are airfoil, with XY. < 1.0, one obtains, from Eqs. 
(18), (20), and (29); 


= Me? — 1)(a2 — a) + 

— &) (34) 
The appropriate change of the angle of attack is ob- 
tained from the equality C, = C, and Eqs. (17), (19), 
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FUEL CONSUMPTION RATIO 7 
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Fic. 9. Effect of body drag on fuel consumption ratio. 


and (29). This gives 


= a — e(ECyW Mo? — (35) 


So that, upon substitution and simplification, one fi- 
nally obtains 


= + a) (36) 


Since ¢,, is of O(1), Eq. (36) indicates that «, is much 
smaller than unity. Consequently, this mode of utiliz- 
ing the wing heat addition should result in energy re- 
quirements which are greater than for a comparable 
airplane without the external heat addition. The dis- 
advantage so incurred will become somewhat smaller at 
high angles of attack and for thick airfoils. 

Therefore, as has been shown, one should not expect 
any substantial reduction in energy requirements re- 


Forces Exerted on a Rigid Wing 


traveling (1) subsonically into the shock, (2) subsoni- 
cally away from the shock, and (3) supersonically into 
the shock. The solution was also obtained for the 
shock striking a stationary wedge. Since experimental 
data are largely absent, the formulas provide a practical 
means of estimating the forces and moments on a wing 
produced by a blast wave. 

The results are also being extended to the case where 
the shock wave is not parallel to the edge of the half 
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by a Shock Wave (Continued from page 80) 


Panel Flutter in the Presence of a Boundary Layer (Continued from page 93) 


sulting from the application of the wing heat addition. 
This was demonstrated when either the wing area or 
the angle of attack are decreased, so as to take the ad- 
vantage of the additional pressure generated by the 
heat addition. Out of these two methods, the reduc- 
tion of the wing area at same angle of attack—e.g., 
redesign of the airplane—proved considerably more 
profitable from the energy requirements standpoint. 
This is so because, when wing area is decreased, the 
energy requirements are reduced due to both the 
smaller pressure and frictional drag. However, when 
wing area is unchanged, all the gains must come from 
the reduction of pressure drag alone. 
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SUMMARY 


Most practical electrical simulators of structures are based on 
the simulation of the relationships of force and moment equilib- 
rium by current-continuity relationships, and on the simulation 
of compatibility of deformations by compatibility of voltages. 
Although equilibrium is readily simulated, the identification of 
fi the voltages with analogous deformations ordinarily requires 
relatively ingenious and involved theoretical analysis for each 
new type of structure. 

This type of simulator is compared in several illustrative cases 
with a second type in which the deformations and their com- 
: patibility are not considered explicitly, but instead the simula- 
tion of the energy in the structure is provided for by resistive 
power loss in the static case and by inductive and capacitive 
energy in the dynamic case. The validity of this type of analog 
depends on the direct use of certain electrical energy theorems 
that are similar to structural energy theorems, with the result 
that the simulator will by itself automatically compute and 
satisfy the compatibility requirements. This use is not to be 
confused with the use of energy theorems in the structure alone 
for the purpose of obtaining deformations preparatory to apply- 
ing compatibility considerations. 

In most structures thus far considered the second type of ap- 
proach has led to a drastically simpler theoretical treatment than 
the first type. In addition, since the analog circuitry is in general 
not unique, the energy approach often leads to a circuit requiring 
fewer components than the compatibility approach and having 
less sensitivity to imperfections of the components. 


(1) INTRODUCTION 


A’ ELECTRICAL SIMULATOR of structures may be 
defined as a circuit which consists essentially of 
passive elements such as resistors, capacitors, and 
transformers arranged in small groups, each of which 
directly represents a single identifiable structural 
member, and in which the currents and voltages di- 
ie rectly represent individual structural parameters such 
as force, moment, deflection, slope, and, in the dynamic 
ik case, the time derivative of these quantities. For the 
present purpose the simulator should be distinguished 
from the older and more conventional analog com- 
s puter in which active members—that is, amplifiers— 
are used in conjunction with passive elements in such 
a way as to represent the derived equations of a given 
ty structure rather than the individual members them- 
selves. Thus, while the computer circuit generally 
bears a pictorial similarity to the equations, the simu- 
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lator bears a pictorial similarity to the structure itself, 
a property which is of considerable usefulness in prac- 
tical applications. 

A good number of simulators have been described.!~* 
In particular, reference 2 cites a large number of papers 
which describe the varied uses of a simulator which 
seems to have enjoyed the widest practical application 
thus far. The basic principles of the analogy used in 
this simulator can be stated in general terms as follows: 


(1) The relationships of equilibrium of forces and 
moments (including inertia force and moment in the 
dynamic case) in the structure are simulated by the 
current-continuity relationships (Kirchhoff’s first law) 
in the analog. 

(2) The contributions of each member to deflection 
and slope are simulated by voltage differences. 

(3) The relationships of overall compatibility of 
slopes and deflections are simulated by the overall volt- 
age-compatibility relationships (Kirchhoff’s second 
law) in the analog. 


These three conditions suffice to produce a valid 
analogy between the circuit and the structure. 

An alternative approach to the problem of simula- 
tion may be obtained by making use of the fact that 
the energy in certain electrical networks obeys laws 
similar to those governing the energy in structures. 
For example, it has been shown‘ that in networks of 
resistors and (practically ideal) transformers, ener- 
gized by externally-supplied constant currents all of 
which are at a single frequency and at the same phase 
(except for 180° reversal to indicate negative polarity), 
there applies a principle of least power which may be 
stated thus: the internal currents assume such values 
as to minimize the resistive power loss, subject to the 
requirements of current-continuity. This parallels 
the case of linear static structures acted on by pre- 
scribed external forces and moments, where the second 
of Castigliano’s well-known theorems (the principle of 
least work) dictates that the internal forces and mo- 
ments assume such values as to minimize the elastic 
energy, subject to equilibrium requirements. This has 
been put to practical application in reference 3, where 
a fourth condition is used in place of items (2) and (3) 
above, as follows: 


7 An alternative arrangement in which the equilibrium rela- 
tionships are simulated by the voltage-compatibility relationships 
will not be considered here, because its use results in a loss of pic- 
torial similarity between the analog and the structure. 
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TABLE 1 
Analogous Parameters in Compatibility and Energy Analog 
(Static Case ) 


Compatibility 


Structure Analog Energy Analog 

force F current / current J 

moment M current J current J 

deflection D voltage V voltage V 

slope @ voltage V voltage V 

work, or elastic (not considered) power + 2, P/2 
energy, W 

force compliance Cr = resistance R = V/] resistance R = P/J? 
D/F = 2W/F? 

moment compliance Cy = resistance R = V/J resistance R = P/J? 


6/M = 2W/M? 


TABLE 2 
Comparison of Necessary Considerations in Compatibility and 
Energy Methods (Static Case) 


Compatibility 


Structure Method Energy Method 
(1) Equilibrium of current-continuity current-continuity 
forces and moments 
(2) Deflection and slope voltage difference 
difference across across individual not considered, 
individual members members but nevertheless 
3) Compatibility of compatibility of obtained auto- 
deflections and voltages matically 
slopes 
t) Work, or elastic (not considered) power + 2 
energy 


(4) The elastic energy of the structure, considered as 
a function only of forces, moments, and elastic con- 
stants, is simulated by half the resistive power loss in 
the analog. 


Thus the structural-electrical analogy is established 
without explicit considerations of deflections and slopes. 

It has further been shown‘ that the electrical circuit 
obeys a law similar to Castigliano’s first theorem, which 
may be stated as follows: in a linear static structure 
acted on by prescribed external forces and moments, 
the deflection (or slope) of any point is the derivative 
of the elastic energy with respect to an external force 
(or moment) applied to that point in the direction of 
given deflection (or slope). The corresponding electri- 
cal theorem, as applied to the circuit described above, 
states that the voltage at the point of entry of an ex- 
ternal current, measured above an arbitrary ground, is 
the derivative of half the resistive power loss with re- 
spect to that current. As a result of this similarity, 
the structural deflections and slopes are automatically 
simulated in the analogous electrical circuit by corre- 
sponding voltages, even though deflections, slopes, 
and voltages were not considered explicitly in setting 
up the circuit according to principles (1) and (4) 
above. 

The technique of simulation using principles (1), 
(2), and (3) above may be designated for convenience as 
the compatibility method, and that using (1) and (4) 
as the energy method. A table of analogous parameters 
with appropriate dimensional symbols is given in Table 
1 for the static case, while Table 2 compares the required 
considerations in the two different methods. 


rigid plate 
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Fic. 1. Rudimentary example. (a) Jeft: structure. (b) right: 
analog. 


A Rudimentary Example 


Fig. 1(a) shows a structure consisting of short con- 
centric tubular members | and 2, compressed by a force 
F acting through an infinitely rigid plate. The mem- 
bers are homogeneous, and have length 1, modulus of 
elasticity E, and respective cross-sectional areas A, and 
A;. The analog is shown in Fig. 1(b), and in this case 
serves for both the compatibility and energy method. 
Following the pattern of Table 2, and noting the 
designation of parameters in the sketch, we write the 
steps in the two methods in a form suited for easy com- 
parison below. (For simplicity here and elsewhere in 
this paper, we ignore the scale factors between elec- 
trical and structural quantities.) 


Compatibility 
Structure Method Energy Method 
(1) F=A+F T=h+h l=h+h 
(2) AD, = L/A,\E = (not considered ) 


AD; = Fe L/AsE = [eRe 


(3) AD, = AD, AV: = AV. (not considered ) 


(4) W = FL/2A\E (not considered ) P/2= 
+ + 


It is also of interest to compare the analytical solu- 
tions of both structure and analog, using both com- 
patibility and energy methods. In the former method, 
steps (1), (2), and (3) lead to the values of resistance 
indicated above, and also yield 


= TR, + Ro), Ts = TR, (R; + Rs) 
AV, = = IR,R2/(Ri + R2) 


and similarly for the analogous structural parameters. 
In the energy method, which we shall apply for the 
analog alone since the structural analysis is very similar, 
the values of resistance are immediately obtained from 
step (4). However, the solution for the currents must 
be obtained by minimizing the power with respect to 
either J; or Js, the other unknown current being de- 
fined by the current-continuity equation. Thus, we 
may write 


OP = 0 = (d, + R| 


This leads to the same current values as before. 
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Fic. 2 (/eft): Truss member. Fic. 3. (right): Analog of truss 
member. 


With the currents known, the voltage differences can 
easily be obtained by borrowing step (2) from the 
compatibility method, but an alternative procedure is 
to use the electrical counterpart of Castigliano’s second 


theorem—tkus, 


AV, = = 
(0/07) + (122R2/2)] 


AV; = 


When the values of /; and /, are inserted, the previous 
values of voltage are confirmed. Naturally, in the 
practical use of the analog, the above differentiations, 
which in more complicated cases can be extremely 
tedious, need not be carried out explicitly. In a cer- 
tain sense, therefore, it may be stated that the analog 
performs these differentiations automatically. 

In more complicated cases than the above the two 
methods of devising the analog may exhibit strong 
differences, with the consequence that under given prac- 
tical circumstances one may be strongly preferable to 
the other. So far as is known, however, these differ- 
ences have not been explored systematically, since the 
various papers describing the simulators have dealt 
primarily with the one or the other. It is the aim of 
the present paper to present certain ideas which may 
help to supply this lack. Mest of these ideas will be 
presented in connection with the synthesis of analog 
circuits for different types of structures, since a greater 
ease of synthesis will lead to wider understanding, 
acceptance, and use of simulators. Following this 
will be given a general discussion of practical consider- 
ations. For simplicity the argument will first be 
referred to the static case, after which the application 
to the dynamic case will be discussed briefly. 

A final word before proceeding to the detailed argu- 
ment—Castigliano’s theorems are frequently used as 
an analytical adjunct to the compatibility method 
because they are convenient in furnishing the relation- 
ships between force, moment, deflection, and slope to 
be satisfied by the analog. For example, the values for 
AD, in the above example might well have been de- 
rived with the aid of Castigliano’s first theorem as the 
partial derivative with respect to F; of the work of 
member 1; namely, F\°L/2A,£. This does not con- 
stitute the energy method under the present definition, 
since it does not make use of the fact that the analog 
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itself obeys laws which are the counterparts of Castigli- 
ano’s theorems. 


(2) Prn-CONNECTED TRUSS 


Fig. 2 shows member a~—b of a pin-connected truss, 
in which the force F, supplied by external and con- 
tiguous members, is necessarily axial. For future con- 
venience F is considered positive when compressive, 
and is presumed not large enough to cause appreciable 
column buckling. The axial compliance is C and the 
member is (initially) inclined at angle ¢ from the indi- 
cated x axis. 

An analog circuit, first presented by Bush,! is given 
in Fig. 3, and is shown in two parts, the x portion simu- 
lating the x projection of force and deflection in the 
member, and correspondingly for the y portion. The 
currents, supplied by external current sources and cir- 
cuits corresponding to contiguous members, are desig- 
nated as shown. The transformer and the resistors are 
considered ideal and the transformer turns are such as 
to give the voltage-difference ratio 1:r, with relative 
polarity as indicated by the plus signs. We are to 
find r, and J, such that the analogy is valid. 

We first establish by force resolution that the currents 
I, and J, may be written as J cos ¢g and / sin ¢, respec- 
tively, where / is a fictitious current which simulatcs 
the real force F. Hence, since the ratio of currents in 
an ideal transformer is known to be the negative re- 
ciprocal of the ratio of voltage differences, we have 


r = —I,/I, = —cot ¢ 


Compatibility Method 

To evaluate the resistances by the compatibility 
method we must first calculate (D,,), and (D,»),, repre- 
senting the respective x and y components of deflection 
Dy, accumulated along the member from a to b. Not- 
ing that D,, is caused both by axial deformation and by 
rotation of the member through a very small angle y 
(assumed positive when counterclockwise), we write as 
a close approximation 


(Dw)r = —FC cos ¢ — Ly sin ¢ 
(Dw)y = —FCsin ¢ + Ly cos ¢ 


where L is the length of the member. The correspond- 
ing expressions for voltage rise in the analog are 


(Vw)r = —IR, cos ¢ + FE, | 
(Vary —IR, sin + rE. 
—IR, sin ¢ — E, cot of 


to 


where £, is the voltage rise in the x winding of the trans- 
former. Hence the analogy will be valid if the resist- 
ances are made equal to each other and proportional 
to C, and £, is interpreted as —Ly sin ¢. It then fol- 
lows that the voltage difference across the x resistor 
simulates the x projection of the axial deformation, 
and the voltage difference across the x winding of the 
transformer simulates the x deflection due to rotation 
of the member; and correspondingly for the y circuit. 
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Energy Method 

Since the elastic energy is F*C/2, we make the re- 
sistive power loss P equal to twice this quantity (accord- 
ing to Table 1), or 


F°C = 1°(R, cos? ¢ + R, sin? ¢) 


so that we are immediately led to the previous values 
of resistance—namely, R, = RK, = C. As to the de- 
flections, we cannot apply Castigliano’s first theorem 
with the assumption that F, and F, are independent 
applied forces, since they are related by equilibrium 
considerations, nor can we consider J, and I, to be 
independent supplied currents, for they are related by 
current-continuity considerations. Instead, however, 
we shall assume that the member of Fig. 2 and the 
x and y networks of Fig. 3 are suitably connected with 
other similar members to form a complete truss. Then 
if we apply load components F,, and F,, (which may be 
fictitious) at a, the deflection components D,, and D,,, 
measured with respect to ground, are given by Cas- 
tigliano’s first theorem, and likewise for D,, and Dy. 
Similarly, supplying the load currents J,, and I,,, we 
obtain, by use of the electrical counterpart of Castigli- 
ano’s first theorem, voltages corresponding to the pre- 
vious deflection components. Thus, it follows that 
the components of deflection difference across the 
member, that is, (D,,), aad (D,»),, are simulated by the 
corresponding voltage differences and (Va»)y- 

Since the voltage differences across the resistors, that 
is, /,R, and /,R,, must simulate the deflection com- 
ponents F,C and F,C due to axial deformation, it follows 
that the transformer voltage differences, which represent 
the differences between total deflection and that due 
to axial deformation, must simulate the components of 
deflection due to the rotation of the member. These 
conclusions have been reached without explicit con- 
sideration of LW and E,, which are of no direct interest 
in themselves. 

By energy considerations, one of the resistors, say 
R,, can be omitted, providing the other is increased to 
a value R,’ such that the total resistive power loss is 
That is, 


I2R,’ = PC 
R,’ = C/cos* ¢ 


unchanged. 


(3) 


The corresponding proof by compatibility consider- 
ations is somewhat more difficult, since it must deal 
explicitly with the altered transformer voltage rise 
E,’. Making use of Eqs. (2), we write 


+ E, = —I,R,' + E,’ 
—1,C + rE, = 0+ rE,’ 


(Vao)z 

(Van)y 
from which £,’ can be eliminated and the preceding 
result for R,’ confirmed. 


Complete Truss 


An example is given in Fig. 4, in which (a) shows a 
simple pin-connected truss loaded by force F, inclined 
at a given direction from the indicated x and y coordi- 
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nates. Figs. 4(b) and (c) show the x and y networks 
of the analog, with the supplied currents designated as 
the components F, and F, of F, and with certain of the 
points grounded (voltage held to zero) to represent the 
fact that the corresponding deflection components are 
held to zero by the truss supports. The branches fol- 
low the pattern of Fig. 3, with the transformer windings 
represented by small circles for convenience, and with 
the electromagnetic coupling indicated by dotted 
lines. Branches 1—2 and 3—4 are omitted in the x net- 
work because they have no x projection. For the sake 
of economy, the resistors in branches 2-3 and 2—4 of 
the x network have beea omitted and the corresponding 
x resistors adjusted according to Eq. (3). The inter- 
ested reader may work out and compare the analytical 
solutions for the structure and the analog, both by the 
compatibility method and the energy method. 


(3) RIGIDLY CONNECTED FRAME 


Fig. 5 shows a member of a rigidly connected frame, 
with end-shears and moments applied by external loads 
and by contiguous members. The reference direction 
for applied moment is counterclockwise, and that for 
shear is such as to tend to rotate the member counter- 
clockwise. Purely for simplicity we shall assume that 
the member is initially straight, is unloaded except at 
the end points (or else the values of shear at the two 
ends might not be equal), and has a constant value of 
I, the moment of inertia of the cross section about the 
neutral axis in bending. The deformation, and conse- 
quently the work, associated with direct shear (and 


Truss and analog. (a) top: truss. (b) bottom left: 
x network. (c) bottom right: y network. 
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Fic. 6 (bottom). 
frame member. 


Fic. 5 (top). Frame member. 


with axial force, if it exists) is considered negligible 
compared to the corresponding bending effects, a con- 
dition often realized in practice. 

A widely used analog due to Russell’ is shown in 
Fig. 6, where the currents /,, M,, and S, representing 
the corresponding structural parameters, are supplied 
by external current sources and by circuits correspond- 
ing to contiguous members. For convenience the two 
separate portions of the circuit may be called the 
moment network and the shear network, respectively. 
The voltage differences on the two sides of each (ideal) 
transformer are designated in terms of the voltage 
differences /,, Es, and the turns ratios 7; and 72, with 
plus signs to indicate relative polarity. The associated 
currents are written in terms of S, 7;, and 7. We are 
to determine the values of 7; and 72 and of the resistors 
R, and R, such that the analogy is valid. 

By equilibrium of moments in Fig. 5, we have 


M,+ M,+ SL =0 (4) 


where L is the length of the member. From this it 
follows by current-continuity considerations in Fig. 6 
that 


(1/n) + (l/r) = L (5) 


Compatibility Method 


We now assume that the member is temporarily fixed 
at initial values of deflection and slope at point a, 
and first evaluate 6,,, which is the increase in elastic 
slope (measured counterclockwise) between a and b, 
and D’,,, which is the increase in deflection between 
a and b, measured from the extended tangent through 
a in the direction of positive shear at b. One conven- 
tional method of determining these quantities is by use 
of Castigliano’s first theorem. (As explained in the 
Introduction, this does not constitute the energy method 
as defined here.) We readily derive the elastic energy, 
or work, in terms of 1, and S, as 


W = (L/2E1) [M,2 + M,SL + (SL)2/3] 


where F is the modulus of elasticity. Differentiating 
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this with respect to \/, and S in turn, we get 
6 = (L/EI) (M, + SL/2) (7) 
Dw = (L?/EI) (M,/2 + SL/3) (8) 
The total deflection change from a to b in the stated 
direction, D,», is greater than the above quantity by 
the amount 6,1, where 0, is the slope ata. Hence 
Dy = OL + (L2/EI) (M,/2 + SL/3) 


In the analog we temporarily assume that points a 
and c are fixed at initial values of voltage, and deter- 
mine the voltage rises V,, and V4 as 


Ves + 5) ro) (10) 

Vea = E, + RS +- E» (1 1) 

Comparing Eqs. (10) and (7), and making use of Eq. 
(5), we find 

R, = L/EI, 

As for Eq. (11), we first note that when Fig. 6 is 

suitably connected with other similar circuits to form 

the analog of a complete rigidly connected frame, the 

voltage at a, 7,4, should simulate the slope 6,, so that 

the voltage at b, 7.F2, will then correctly simulate 4). 

Then by the transformer arrangement of Fig. 6 we have 


= 1/r = L/2 


Via RS + (L 2) (0, 6y) = 
RS + (L/2) (20. + 40) 


which with the aid of Eq. (7) can be seen to simulate 
Eq. (9) if 


R, = L*/12EI 
thus completing the analogy. 


Energy Method 


The expressions for 7; and 72 are obtained in the same 
manner as before. To evaluate the resistances we 
write the resistive power loss as 


P = RS? + Ri(M, + SL/2)?* 


On equating this to twice the value in Eq. (6) we 
quickly obtain the same values of resistance as before. 
Further, by Castigliano’s first theorem and its electrical 
counterpart, when Fig. 6 is suitably connected to other 
similar circuits to form the analog of a complete 
rigidly connected frame, the slopes at a and b in the 
structure will be simulated by the respective voltages 
at a and b in the analog, because external moments 
applied at a and b in the structure would be simulated 
by external currents supplied to a and b, respectively, in 
the analog; and deflection components at a and b in the 
structure will be simulated by the voltage at respective 
points in the analog where external currents would be 
supplied to simulate corresponding components of 
external force in the structure. Thus, all quantities of 


major interest can be simulated without considering 
slopes and deflections explicitly in the process of setting 
up the analog. 
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Complete Frame 


As a simple example, consider the rigidly connected 
frame of Fig. 7, consisting of two members as in Fig. 
5. The applied load is F, and the force and moment 
reactions are as indicated (the D’s should be ignored 
for the moment). For convenience the shear in each 
member is designated. 

One form of the analog, shown in Fig. 8, can be de- 
vised by first connecting together the individual mo- 
ment networks at joint 2, grounding the combined 
network at 1 to simulate zero slope and opening it at 
3 to simulate zero moment. As before, each trans- 
former winding is symbolized by a small circle as shown 
in the Figure, the two circles corresponding to each 
transformer being connected by a dotted line to indi- 
cate the electromagnetic coupling. The individual 
shear networks may then be drawn, and are first left 
open at all ends except at joint 1’, which is grounded to 
simulate zero deflection. (Certain lines have been 
shown dotted so that the portion of the analog thus far 
described can be visualized easily.) The currents 
traversing all transformer windings thus far drawn are 
indicated, thereby fixing the transformer turns ratios, 
which are not shown explicitly. Considering also the 
current /;, we may now deduce by overall current- 
continuity considerations that 


FLy = 0 


Rigidly connected frame. Fic. 8 (bottom). Com- 


Fic. 7 (top). 
patibility analog of frame. 
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Fic. 9. 


Energy analog of frame. 


which is precisely the equation of overall moment- 
equilibrium. 

Next we shall connect the individual shear-networks 
according to compatibility considerations. Desig- 
nating D>_,2, for example, as the deflection, with respect 
to ground, of joint 2 (indicated by the first subscript) in 
the direction of positive shear at end 2 (indicated by 
the last subscript) of member 1-2, and similarly for 
certain other deflections, it is possible to deduce by 
purely geometrical considerations that 


= 
= Ds_12 


These various quantities are indicated as appropriate 
voltages at points 2’ 2”, and 3’ of Fig. 8. The above 
relationships are enforced by the dotted lines associated 
with transformers 7; and 7», whose turns ratios are 
indicated, with relative polarity as given by plus signs. 
As a result of these turns ratios the shear current 
flowing in 2”-3’ induces the two currents H 2, which 
together with the current —F supplied by an external 
source as shown satisfies the current requirement in 
1’-2’.. The analog is now complete. 

Next we develop the analog by energy considerations. 
This will be discussed in connection with Fig. 9, where 
some lines are dotted so that the steps used in the con- 
struction of the circuit may be easily visualized. 

First we connect individual moment networks as 
before, ignoring the transformers and showing only the 
currents and resistors, combining the two transformer 
currents entering joint 2 to form —FL,/2; this leaves 
unchanged the similarity between the equations of 
moment-equilibrium and the corresponding equations 
of current-continuity. Of the shear network for mem- 
ber 1-2 we show only the resistance, and multiply the 
current by the factor L)/2 so that, insofar as current- 
continuity is concerned, point a can be directly con- 
nected to a’. Similarly, in the shear network for 
member 2-3 the current is multiplied by the factor 
—I/+/2 so that, insofar as current-continuity is con- 
cerned, points b and b’ can be directly connected. 
The two shear resistances must be multiplied by the 
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Rectangular panel. Fic. 11 (bottom). Analog of 


Fic. 10 (top). 
rectangular panel. 

respective factors 4/L,? and 2/Z,? so as to keep the 

resistive power loss unchanged. 

Next we must merely relate the currents traversing 
points c, d, and e in a consistent manner. This can be 
done as shown by the dotted lines, using a current 
source to introduce — FL, and a 1:1 transformer, with 
relative turns polarity as indicated by the plus signs, 
to divide this current into two equal parts. The analogy 
is now complete. 

By Castigliano’s first theorem the joint slopes are 
simulated by joint voltage as before, while the voltage 
at the point of entry of — Fly) represents — (—Dz2_12)/Lo 
as shown. D3 23; is not simulated directly, but the 
quantity D;,, the deflection of point 3 in the vertically 
upward direction (which is actually of more direct 
interest than D3 »;), can be simulated, if desired, by 
imagining an upward load V applied at point 3, re- 
arranging the analog to simulate this load, and noting 
the voltage at the point of entry where current would 
be introduced if the load were actually present. With 
the aid of the same principles as before this modifica- 
tion of the analog could readily be accomplished by 
means of a source of current VZp and a 1:1 transformer, 
as shown in the Figure, connecting point d’ to d and 
e’ toe. In fact, merely by inspection of transformers 
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T, and 7. we deduce with the aid of Castigliano’s theo- 
rem that 
D3, Lo D» 12 Ly 


which could readily have been deduced by purely 
geometrical considerations. 

It may seem surprising that the analogs in Figs. § 
and 9, one of which is so much simpler than the other, 
should furnish the same values of the significant param- 
eters. That this is so, however, can readily be checked 
in the particular case where 


= The, 


Exact values of some of the resulting parameters are 
as follows: 


H = — F/6, 
(3, 5) Oe, 


M, = = (5 12) FL? Ely, 
Ds_12 = —(11 36) FL Ely 


With these values the unknown transformer voltage 
drops (which are of no direct interest in themselves) 
may be readily evaluated, taking due account of the 
magnitude and polarity of turns ratios, to permit the 
detailed check of Figs. 8 and 9. 


(4) STIFFENED FLAT SHEET 


Fig. 10 shows a typical internal rectangular panel j, i 
of a stiffened flat sheet, loaded at external boundaries 
only, of a general type for which electric analogs have 
been discussed by Newton,® Ross,’ and Goran.* The 
y directed members, or stringers, have cross-sectional 
areas A which are assumed large enough compared to 
the nearby cross-sectional area of the sheet so that these 
members resist virtually all the force in the y direction, 
and are further assumed to have negligible extent in the 
x direction. The x directed members, or ribs, are 
shown dotted because they may be fictitious, and are 
presumed to represent the resistance of the sheet to 
x force; these members are presumed to be infinitely 
stiff and to have negligible extent in the y direction. 
The panel sheet is presumed to have a constant thick- 
ness ¢;, ,and to be subject to constant shear stress acting 
in the directions of the stringers and ribs. Buckling is 
presumed negligible. Though rather coarse, these 
assumptions have often been considered useful, as for 
example by Duberg.’ 

The tensions in the stringers and ribs are shown as P 
and Q, respectively, the lengths as L, and the shear 
flows in the panel and in adjacent panels, considered 
in the direction in which they exert force on the ribs 
and stringers, are shown as g. Considering the equi- 
librium of y and x forces, we have 


— Py = Li 


i) = 0 (12) 


gi, i-1) = O 


Newton and Ross do not deal explicitly with the 
second of these equations, with the consequence that 
the application of their analog is limited to those par- 
ticular cases in which this equation does not influence 
the solutions, as is pointed out on p. 411 of reference 5. 
The latter overcomes this limitation by writing the 
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equations of moment-equilibrium, which when used 
in conjunction with the y force equilibrium equations 
generally supply the same information as the x force 
equilibrium equations. 

To obtain the compatibility relationships, Goran ex- 
presses the elastic energy of the panels, consisting of 
tensile energy in the stringers and shear energy in the 
sheets, and minimizes this analytically, subject to the 
equations of equilibrium of y forces and moments. He 
thus obtains certain relationships which he then satis- 
fies by voltage-compatibility considerations in the ana- 
log. By the present definition, this constitutes a com- 
patibility method. The resulting analog is shown in 
Fig. 11, in which the equations of current-continuity 
simulate the first of Eqs. (12), and the small circles 
represent the primaries of transformers having turns 
ratios inversely proportional to the moment arms of the 
forces P;,,,; about some suitably designated point, 
so that the secondaries can be connected in such a way 
that the resulting equations of current-continuity 
simulate those of moment-equilibrium. On the as- 
sumption that P;, ; and A,;, ; remain constant between 
points j, 7 and j, 7 + 1, the resistors are found to have 
the following values: 


R,, = L,/A;, E 
L;/Lit;, G 


where E and G are the moduli of elasticity in tension 
and shear, respectively. 

Proceeding by the energy method, 
equilibrium relationships are satisfied, whether by 
Goran's network or any other, it is merely necessary to 
choose the resistance values so that the electrical power 
equals twice the elastic energy. By this approach one 
may easily prove the validity of Goran’s analog without 
following his rather elaborate, though interesting, 
compatibility analysis. 

The interested reader may refer to Goran for further 
details as to the connections of the secondaries of the 
transformers in Fig. 11, and for means of simulating 
the boundary conditions in problems involving com- 
plete stiffened plates. 

The advantages of easy circuit synthesis obtainable 
by the energy method are likely to be increasingly im- 
portant as the rough idealizations cited at the beginning 
of the present discussion are progressively refined. 
For example, if allowance is made for the fact that the 
transverse stiffeners (or their equivalent, obtained by 
assuming a certain portion of the plate concentrated 
along each stiffener axis) have finite rather than in- 
finite stiffness, this is readily taken into account in the 
analog by inserting a resistance of appropriate value in 
the path of a current which simulates the force in each 
stiffener. 


once the 


(5) PrRacTicat CONSIDERATIONS 


Circuit Economy 


Obviously there are cases where the energy analog 
requires significantly fewer components than the com- 
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patibility analog, as can be seen by compariscn of Figs. 
Sand 9. Such a saving in components is of practical 
importance in the solution of complicated problems, 
where the available number of components is limited by 
cost considerations. While the economy of Fig. 9 
could possibly be achieved by suitable application of 
compatibility considerations to Fig. 8, this would re- 
quire relatively difficult and ingenious supplementary 
analysis because the deflections and slopes would have 
to be considered explicitly. 


Analog Error 


Since the components are never ideal, the decrease 
in the number of components which is sometimes 
yielded by energy considerations may generally be ex- 
pected to lead to decreased error. However, there 
are more important ways in which energy considerations 
can be used to decrease error, as will now be discussed. 

As a practical matter the principal sources of error 
are the series resistances of the transformer windings 
and the equivalent transformer shunt resistances, each 
of the latter being referred to either winding of its as- 
sociated transformer. (The equivalent series and 
shunt reactances of the transformers, which may well 
be of the same order of magnitude as the winding and 
shunt resistances respectively, will not introduce signif- 
icant error if the analog currents and voltages are 
measured by certain well-known devices which are 
sensitive only to the in-phase components of these 
quantities.) In both the compatibility and energy 
methods each winding resistance may merely be ab- 
sorbed within a true circuit resistance with which it 
happens to be connected in series. However, not all 
windings are thus connected to true resistors. Never- 
theless, by energy considerations cach winding re- 
sistance can still be compensated perfectly if a resistor 
is found anywhere tn the analog carrying a current pro- 
portional to that in the winding, this resistor merely 
being decreased in value by such an amount as to keep 
the total power loss unchanged. To determine and 
justify such a compensation by compatibility consider- 
ations is often impractically tedious. 

As to the equivalent shunt resistances, these cannot 
ordinarily be compensated in advance. However, the 
voltages to which such shunt resistances are exposed 
may differ radically in the compatibility and energy 
analogs. In Fig. 6, for example, the shunt resistances 
of the transformers may be referred in each case to 
the winding connected to the moment network, so 
that the voltage across a typical shunt resistance is 
proportional to the slope at the joint to which it is 
connected; this is seen in Fig. 8. In the energy 
method, however, it turns out to be rather typical for 
each moment network to be fed at its two ends by an 
arrangement similar to that of moment network 1-2 
and transformer 7> in Fig. 9, with the result that, ignor- 
ing the voltage drop in the shear resistance for the pres- 
ent purpose, the voltage across one winding of a typical 
transformer is proportional to about half of the differ- 
ence between the slopes at the two ends of the member. 
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In complicated frames this quantity is on the average 
very much less than the slope at each joint, thus gener- 


: ally leading to smaller error due to transformer shunt 


resistances. 


(6) Dynamic CASE 


In the dynamic case the compatibility method follows 
much the same pattern as before, except that inertia 
forces are taken into account; the resulting analog 
need not be discussed in detail here, except to say that 
it generally employs inductors instead of the previous 
resistors, and uses capacitors to represent masses as- 
sumed concentrated at certain points in the structure, 
each capacitor being connected between ground and the 
point in the analog whose voltage represents the fime 
derivative of deflection of the point at which the corre- 
sponding mass is assumed concentrated. 

The energy method is based on the fact that La- 
grange’s well-known dynamic equations, which in the 
structure are derivable from a consideration of elastic 
and kinetic energy, also apply in electrical circuits of 
inductors, capacitors, and ideal transformers, with 
inductive and capacitive energy taking the respective 
places of elastic and kinetic energy in the structure. 
(Though long known, this fact has recently been proved 
by conventional electrical engineering techniques”) 
With the aid of these ideas it has been shown in a very 
general way’! that any analog of a static structure can 
be changed to a dynamic analog by replacing the re- 
sistors by inductors of proportional value, and by con- 
necting a capacitor between ground and each point 
whose voltage represents a time derivative of deflection, 
the capacitor values being such that the kinetic energy 
is appropriately simulated. 

Thus, the relative advantages and disadvantages 
of the compatibility and energy methods are carried 
over essentially unchanged from the static to the dy- 
namic case, with one exception—in the energy analog 
one must explicitly consider those voltages which rep- 
resent time derivatives of deflection of points where 
mass is assumed concentrated, which may mean a 
partial loss of the circuit economy gained as against 
the compatibility analog in the static case. For ex- 
ample, if in Fig. 7 a mass were assumed concentrated 
at joint 3, it would be necessary to use transformer 7) 
in Fig. 9 in order to develop the voltage D;,/L». 


(7) CONCLUSION 


The foregoing comparison of the compatibility and 
energy methods, the latter implying the explicit use of 
energy principles both in the structure and in the ana- 
log, is by no means exhaustive. Nevertheless, it has 
been shown that the energy method, in which there is 
often no need for explicit consideration of the structural 
deformations and the electrical voltages (though the 


former are usually readily obtainable if required), fre- 
quently has some or all of the advantages listed below. 

(1) The synthesis of the analog circuit is easier since 
fewer parameters are considered. This is important 
in two ways: first, in simplifying the treatment of 
types of structures already considered by the com- 
patibility method; and second, in the extension of simu- 
lation theory to new and difficult cases, as for example 
when the coarse simplifying assumptions cited above in 
connection with stiffened sheets are increasingly refined 
for the sake of greater realism. 

(2) Fewer components are required, thus decreasing 
the cost and complexity of the simulator required to 
handle a given problem. 

(3) The possibilities for compensation (in advance) 
of transformer imperfections are greater, and the cir- 
cuit can be arranged so that uncompensated imper- 
fections have a less adverse effect. Thus, greater 
transformer imperfections may be tolerated, with a 
consequent advantage in cost. 

The advantage of easier circuit synthesis has been 
found to apply in all of the nontrivial cases thus far 
considered by the writer, and in no case has the cir- 
cuit of the energy analog been found more complex 
or more sensitive to component imperfections than 
that of the compatibility analog. Hence, the present 
emphasis on compatibility concepts (including the 
case where energy considerations are used only in the 
structure to obtain deformations) should be critically 
examined, with a view to greater use of the energy 
method. 
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On “‘A Critical Strain Approach to Creep 
Buckling of Plates and Shells’’ 


Nicholas J. Hof 


Head, Division of Aeronautical Engineering, Stanford University, 
Stanford, Calif. 


August 25, 1958 


I* A RECENT PAPER,’ good agreement was shown between cal- 
culated and measured strains at the time of creep buckling, 
and between calculated and measured critical times, in the case of 
circular cylindrical specimens subjected either to pure compression 
or to pure torsion. The conjecture underlying the calculations 
was that an initially perfect specimen would buckle after the lapse 
of sufficient time to permit the development of total creep strains 
equal to the strain calculated from the relevant elastic buckling 
formula from which Young’s modulus was stricken. For instance, 
at the moment of buckling the total creep strain of a uniformly 
compressed circular cylindrical shell was taken as 0.64/R, where 
h is the wall thickness and R the mean radius,' while the total 
creep strain of a simply supported column was given as r?/(L/p)?, 
where L is the length of the column and p the radius of gyration 
of its cross section.” 

As the range in which the parameters of the specimens were 
varied was very restricted, it would be unjustified to assume that 
the article proved the conjecture to be true. Asa matter of fact, 
the writer of the present note was able to prove’ that the assump- 
tion of initial perfection leads to conclusions entirely different 
from those proposed by Dr. Gerard. Moreover, the engineering 
significance of any results derived for initially perfect columns is 
rather doubtful on the basis of reference 3. 

Dr. Gerard’s conjecture appears particularly unjustifiable if one 
realizes that elastic-plastic columns that do not creep also 
do not buckle at the strain values quoted. This is clear if the 
Ramberg-Osgood‘ representation of the stress strain curve is 
used: 
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= (o/E) + (1) 
Jt was shown in a recent paper® that the Euler strain 

eg = (2) 

of a pin-ended column made of such a material is given by 
= (o/E) +1 (o/B)" (3) 
according to the tangent modulus theory when buckling occurs 

during a monotonically increasing loading process, and by 
ég = + (1/2)n (0/B)” (4) 


according to the reduced modulus theory when the load remains 
constant during deformations, as in the creep buckling process. 

Let us assume that the buckling stress of a particular column is 
such that the plastic strain (o/)" is exactly equal to the elastic 
strain (o0/E) in Eq. (1). The value of m for the particular mate- 
rial should be 10. Then the Euler strain of Eq. (3) is five and a half 
times, and that of Eq. (4) is three times the value of the actual 
strain at buckling. But Dr. Gerard suggests that the strain at 
the moment of buckling should reach these larger values. 

As far as the tests described in reference 1 are concerned, it 
should be noted that they were carried out with cylindrical shells 
having radius to thickness ratios of 9.7, 14.1, and 18.4. This i 
not only a very restricted range but also one that is not commonly 
encountered in aeronautical engineering. Moreover, the defini 
tion of the buckling strain as the one corresponding to the inter- 
section point of two straight lines passing through the experimental 
points is rather arbitrary; in the test results known to the writer, 
such a definition would lead to serious difficulties. Finally, the 
selection of 0.6 as the coefficient in the shell buckling formula is 
also debatable. 

The writer agrees with Dr. Gerard that the creep buckling 
problem is still far from reaching a satisfactory solution. It is 
hoped that much more work will be done in this field both theo- 
retically and experimentally. 
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Comments on ‘‘The Determination of the 
Shear Center for a Special Solid 
Symmetrical Airfoil’: 


J. Lockwood Taylor 
Professor, Norges Tekniske Hogskole, Trondheim, Norway 
July 28, 1958 


it Mr. SISTINO’s NOTE,! I find an error in the first term of 
Eq. (5) (right-hand side) in which y should be replaced by 
(y — C), where C is the coordinate of the centroid of area of the 
section, 3/7 in the present instance. Without this correction, 
his solution is for a beam twisted, as well as bent, and the center 
he calculates is, therefore, not the true shear center. Timo- 
shenko himself makes a similar mistake in one of his examples, 
that of the semicircular section. 

It is a mathematical oddity that exact simple solutions of 
the type given by Mr. Sistino are available only for certain 
combinations of thickness, h, and Poisson’s ratio », given the 
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value of C. I find that such a solution can be obtained for 
h = 0.2—i.e., practically Mr. Sistino’s figure—and » = 0.321, 
This gives the shear center at 0.3865 chord from the leading 
edge, against the value he gives of 0.561 chord. It is not to be 
expected that the assumed value of v will be at all critical. 
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Author’s Reply 


Alphonse J. Sistino 

Analytical Engineer, Fairchild Engine and Aircraft Corporation, 
Deer Park, N.Y 

August 28, 1958 


] WIsH to thank Prof. J. Lockwood Taylor for indicating the 
incorrect numerical value of the shear center, as given in 
my paper, ‘Determination of the Shear Center for a Special 
Solid Symmetrical Airfoil.’”, However, I do not agree with his 
suggestion that a correct result can be obtained by replacing 
y by (y — c) in the first term of Eq. (5) (right-hand side) where 
c is the coordinate of the centroid of area of the section. In 
order to arrive at a correct result, the origin of the coordinate 
system should be located at the centroid of the cross section. 
If this is done, then 


x= 
¢ = — + — y — + 
and Eq. (5) is correct. 
Following the mathematical procedure outlined in the paper, 
with the above modifications we have then 
k= 0316 
Ss 0.125 P/I 
r = 0.002 P/I = 0.002 P/I 


and the shear center distauce is 0.375 from the leading edge. 
This result is for an airfoil which is 1.564 thicker than shown in 
Fig. 1 of the original paper. 

Also, I wish to correct a misprint in the original paper. The 
coordinates (x, v) shown in Fig. 1 should be reversed. 


A Rapid Method for Estimating Grit Size for 
Fixing Transition on Wind-Tunnel Models 


Jerold M. Bidwell 
Design Engineer, The Martin Co., Baltimore, Md. 
August 25, 1958 


A CHART HAS BEEN PREPARED (Fig. 1), using the results of two 
well-known papers,’ ? that permits a rapid estimation of the 
grit size required for fixing transition on wind-tunnel models. 
This chart has been prepared for wings. A similar chart can be 
prepared for axisymmetric bodies by using the Mangler trans- 
formation on the Chapman-Rubesin equations.” 

Braslow! and his associates have shown that if the top of the 
grit is within the “‘linear’’ portion of the boundary layer and if 
the roughness Reynolds Number, Ry = (uxy,/ve), is equal to or 
greater than 600, the transition to turbulence will take place at 
the grit element. Braslow also observed that stabilizing the 
boundary layer by pressure gradient or cooling does not affect 
the criterion. 

The Chapman-Rubesin equations for the boundary layer on 
an insulated flat plate are, from reference 2, 
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T/T, = 1+ — 1)/2]M,,? r(n) (36 mod. ) 


(y/2) VU = 0 + — (42 mod.) 


u/U.. = (1/2)f'(») (13) 
= C(T/T,) (4) 


The pressure is constant through the boundary layer, so 

p/p tif, 
By substituting these equations in the Braslow criterion, the 
following two equations evolve, which must be solved simul- 


taneously : 
ye = [12000/(U/v)] + — 19/2] 
ye = [2V xc/(U,/v~)] + — 1)/2] M2 


where the subscript & refers to the values at the top of the grit. 
The solution to these equations is presented in Fig. 1, for a = 

1.4 and C = 0.90. Deviations of about 20 per cent in C give 

deviations in y;,/x of about 10 per cent. The deviations in this 


range may be estimated by the approximate formula 
= + 0.885) 


In the figure, x is the distance from the leading edge to the rough- 
ness strip. The Reynolds Number, R,, is based on x and the 
conditions outside of the boundary layer. For thin wings, this 
is assumed to be free-stream conditions. 

The chart has been prepared primarily for wind tunnels oper- 
ating at stagnation temperatures of about 150°F. Extrapolating 
the curves should be done with caution. For other applications, 
it is recommended that the paper by Braslow be studied care- 
fully.§ 
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A Similarity Law for the Normal Force on a 
Delta Wing at Hypersonic Speeds 


Arthur F. Messiter 

Member of the Technical Staff, Hughes Aircraft Co., 
Culver City, Calif. 

August 25, 1958 


= A FLAT DELTA WING at an angle of attack a in a free 
stream of high Mach Number V/, such that A/ sin a >1. 
If the ratio of specific heats y is assumed close to one, it is ex- 
pected that a reasonable estimate of the normal-force coefficient 
Cy is given by the Newtonian value Cy = 2 sin? a. Using this 
value as a first approximation, one can derive equations describing 
the consistent second approximation, and a similarity law can be 
obtained for Cy. 

In order to allow large angles of attack, @ is chosen to be O(1). 
Also we will take (y — 1) = O(1/M?*). For sufficiently large as- 
pect ratios, the strong shock wave on the windward side is 
attached along the wing leading edge; the distance from the wing 
surface to the shock can be shown to be O(1/.M?). Defining } as 
the tangent of the wing semiapex angle, one can show that, as } 
decreases in order of magnitude, the shock wave detaches from the 
leading edge when b = O(1//). In order that the theory allow 
the shock to be either attached or detached, we must choose ) = 
O(1/M). Thus, the orders of magnitude to be considered are 


M>1, a=O(1), y —1 = O(1/M?), b = O(1/M) (1) 


and 
Mb held fixed. The case ©,a—~0, can be shown 


That is, we consider the limit JJ ~ © with a, (4 


to be a limiting case of Eq. (1). 

The necessary equations include the usual differential equations 
for inviscid flow, the three-dimensional shock relations, and the 
boundary condition specifying no flow through the wing surface. 
The requirement that these equations be nontrivial in the limit 
determines the orders of magnitude of all the variables. Know- 
ing these orders of magnitude, one can then introduce an appro- 
priate set of new variables, each of which is of order unity. Cer- 
tain constants of proportionality are included in the definitions 
of these quantities in order to make full use of the similarity 
properties of the equations. 

The wing vertex is taken as the origin of coordinates. The 
x-coordinate is positive downstream along the wing centerline; 
z is normal to the wing surface and positive on the compression 
side; and the y-axis is perpendicular to the x-z plane, directed so 
that the coordinate system is a right-handed one. In defining the 
distorted coordinates, the conical property of the flow is utilized: 
y* = y/bx, z* = 2/b*x tana (2) 
If the velocity vector behind the shock, made nondimensional 
with respect to the free-stream velocity, is denoted by g/U = 
(cos a + u, v, w), the components uw, v, and w are found to be 
small, and the new perturbation velocities can be defined as 
follows: 


u* = u/b? cosa, v* = v/bcosa, w* = w sin a/b? cos? a (3) 


The pressure and density, referred to free-stream quantities, are 


p/p. U? = sin? a + (1/M?*) + b? cos? a p*(y*, s*) + 
O(1/M*) (4) 


= — 1)/2] + (1/M? sin? a) + O(1/M"*) 


When substituting into the full equations, it is convenient to 
introduce a new parameter m. 


ol 26 28 1.0 x 10 
Ry 
| 


3 
4 


1/m? = (tan? a/b?) { [(y — 1)/2)] + (1/M? sin? a)} (5) 


The quantity m can be interpreted physically as the ratio b/tan y, 
where yu is the Mach angle in the flow behind the shock. 

Neglecting higher order terms, the differential equations for 
the second approximation are as follows: 


+ = 0 


(v* — y*) + (w* — = (6) 
— y*)v*y* + (w* — v*,* = 0 
(v* — y*) w*y* + (w* — 2*) = —(1/m?) p*+ 


The energy equation has been omitted; it serves only to deter- 
mine the perturbation in density. Jump conditions at the shock 
wave can be formulated in the manner of reference 1. If the 
shock shape is given by 2* = S*(y*), the limiting form of the 
shock relations becomes 


u* = —S* + y*(dS*/dy*), v* = —(dS*/dy*) (7) 

w* = —u* — v*? — (1/m?), p* = w* — u* ' 
where all quantities are evaluated at s* = S*, At the wing sur- 
face the boundary condition requires 


w* =Oatz* =0, -1<y*<1 (8) 


The obvious and striking conclusion to be drawn from Eqs. 
(6)-(8) is that the variables u*, v*, w*, and p* depend only on 
the coordinates and on the single similarity parameter m. That 
is, the similarity properties of the equations are such that the four 
parameters M, a, y, and 6 enter the problem only in one particu- 
lar combination, and the problem is really a one-parameter 
problem. 

The meaning of m can be made a little more clear by consider- 
ation of shock detachment. If the shock is attached along the 
leading edge, there is a region of uniform flow near the leading 
edge, where w = 0 and bv = —u. It follows from the shock rela- 
tions Eq. (7)-that 


b/p..U? = sin? a + (1/M?) + (m?/2) sin? {[(y — 1)/2] + 
(1/M? sin? [1 — V1 — (4/m*)] (9) 


For m < 2, the solution no longer exists because the square root is 
imaginary, and hence the theory predicts shock detachment at 
m = 2, Atm = 2 the flow direction at the leading edge behind 
the shock makes an angle b/2 with the leading edge, and the 
velocity component normal to the leading edge is sonic. While 
the conventional meaning of the expression ‘“‘sonic edge’’ does 
not apply directly, one might redefine this expression in terms of 
the velocity behind the shock. In this sense, it could be said 
that the shock detaches when the leading edge is sonic. 

Following the usual practice, the pressure on the lee side is 
set equal to zero. Integration of the pressure over the wing 
then gives 


Cy = 2 sin? a + (2/M*) + sin? af [(y — 1)/2] + 
(1/M? sin? a)} F(m) (10) 


1 
where F(m) = 2 mf p*(y*, 0; m)dy* 


Eq. (10) is the similarity law for Cy. While the first term seems 
generally to be called the Newtonian value, a pressure coefficient 
of 2 sin? a would give the second term as part of the Newtonian 
result. In the present context, 2 sin? a represents the first ap- 
proximation, and the remaining terms are correction terms de- 
rived for the second approximation. 

Experimental data can be correlated according to Eq. (10) in 
order to study the range of validity of the theory. This is cur- 
tently being done, and the correlation is observed to be good for 
a much wider range of M and a@ than would be anticipated— 
namely, for all values such that M sina>1. The form of F(m) 
has been determined from the available data, and it is seen that 
F(m) varies from —1 at m = 0 to +2 at the point of shock at- 
tachment. For larger m, F(m) remains constant and equal to 2, 
and the normal force agrees with the result of two-dimensional 
shock-expansion theory, taken in the Jimit as M— o and y > 1. 
The correlation of experimental data and the derivation of 
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Eqs. (6)-(8) will be described in detail in a subsequent paper, 
Also to be discussed are a method of integrating the differential 
equations and a method of taking into account the wing thickness, 
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On the Contribution of Thermal Creep to 
the Skin Friction in Slip Flow 


A. F. Charwat* 
Associate Professor, Department of Enounion, University of 
California at Los Angeles, Los Angeles, Calif. 


August 21, 1958 


IHE BOUNDARY CONDITIONS for rarefied-gas boundary-layer 
flow are derived by studying the interference between a 
Maxwellian emission of particles from the wall and the stream 
above it. This yields a ‘‘slip’’ velocity at the wall which de- 
pends on the gradient of the streamwise velocity and also on a 
second-order term involving the wall temperature distribution 


u, = 1,[k(Qu/dy) + (@T/dx)] 5 (1) 


where k and k! are numerical constants dependent only on the 
nature of the gas and the subscript s refers to quantities evaluated 
at the wall. 

The ‘‘thermal creep” is a known phenomenon, and it has been 
studied in connection with the flow through tubes and the 
“radiometer” problem. The same effect—i.e., a shear stress 
opposing the wall temperature gradient—is to be expected also 
in the case of strongly heated aerodynamic surfaces. 

A qualitative estimate of the thermal strain on the skin friction 
of a flat plate can be observed by extending the analogy between 
it and the Rayleigh problem studied by Schaaf? and Mirels* who 
considered only the first term of the slip velocity. The problem 
is defined by requiring that the streamwise derivatives of the 
velocity and pressure vanish identically; the flow is nonsteady, 
and the boundary conditions are chosen to imply a plate along 
the x-axis in an initially stationary fluid which is impulsively 
started from rest at the time origin with a constant velocity ™. 
If the Prandtl Number (¢) is constant then the streamwise- 
momentum equation can be reduced (for details see the original 
work or a summary in reference 1) to 


0u/Or = (a/Re,) (0%u/dy?) (2) 


and can be solved independently of the energy equation. +r = 
u,t/l, is a dimensionless time parameter; /, is the mean free path; 
and Re, is the Reynolds Number based on the mean free path at 
a reference condition (subscript 1). y is a stream function 
satisfying the equation of continuity. 

However, the boundary condition (1) at y = 0 is inconsistent 
with the assumption defining the problem, except in the special 
case when the term containing the temperature gradient is con- 
stant. This term is 


Q = y) X 
(3) 


For this particular case the wall temperature varies with distance 
along the plate as the 3/2 power: 


(T./T1)3/2 = 16V y M(P./P1) (x/h)Q + K (4) 


where K is an integration constant, the temperature of the plate 
at the origin of the coordinate system. 

Eq. (2) can be integrated by standard methods in the r — y co- 
ordinate system. It is then possible to evaluate quite simply the 
local friction coefficient and the wall slip velocity at y = 0 with- 
out having to transform the solution to the physical plane. The 


* The author would like to acknowledge the help provided by R. B. Babco. 
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results are 
Cp = (20/Re,) (1 — Q) exp (or/Re;) erfe (ar/Rei (5) 
U./u = + (1 — Q) exp erfe (a7/Re,)'/? (6) 
when Q = 0—..e., the wall is isothermal— Eqs. (5) and (6) revert 
to the solution of reference 3. 

The ‘‘analogy’’ between this nonstationary problem and the 
semi-infinite flat plate consists in setting the time from the start 
of the motion proportional to the distance from the leading edge 
of the stationary plate. In the present problem it is necessary 
to simultaneously identify the leading edge of the stationary 
plate with a point—for instance, the origin of the nonstationary 
plate. If the first assumption is true, which implies that the 
boundary layer has no ‘“‘memory”’ of its past development, then 
the second one appears to be also acceptable because the Ray- 
leigh flow along the moving plate is identical over all sections at 
any time. The solution is physically meaningful only in the 
region where the temperature is finite and positive. However, 
the ‘‘analogy’’ between the stationary and nonstationary problem 
under the transformation t = Ux is weak. As the inspection of 
the complete set of the governing equations reveals, the condi- 
tions for equivalence (such as v = Q) are not rigorously valid 
even for the isothermal-wall compressible case.* (They are valid 
for the incompressible flow.?) The addition of a streamwise 
variation of wall temperature—consequently a streamwise vari- 
ation in density (the pressure was forced to be constant )—and, 
by continuity, the existence of a cross-flow velocity make the 
analogy even more approximate. 

In spite of the shortcoming of this analysis, a qualitative esti- 
mate may be attempted. There is no reason to foresee that the 
trend indicated by this model should be incorrect. The ratio of 
skin friction with ark without wall temperature gradient is 
Cry/Cr, = (1 — Q) = 1 — 4.24 K 1073 (1/uPg) X 

7,3/?)/dx] 


where Cr, is simply given by Eq. (5) with Q zero. 

The numerical constant was calculated using the ideal-gas 
equation of state and the values of gas constants characteristic 
of diatomic air. The units are in degrees Kelvin, centimeters, 
and atmospheres. This ratio exceeds unity when the temper- 
ature falls downstream of the leading edge. At typical altitudes 
(for instance: P = 10° atm.) and for high temperature gra- 
dients (for instance: 20°K. per centimeter at 1,000°K.) the 
“thermal strain’’ increases the skin friction by about 5 per cent 
at 20,000 ft./sec. flight speed. 

The result is quite sensitive to the assumption of temperature 
gradient and ambient pressure. The thermal stress term can be 
as high as 50 per cent of the slip-flow friction, if a higher temper- 
ature gradient or a lower pressure were assumed. However, if 
one considers an aerodynamic surface, the influence of these fac- 
tors tends to balance out; higher pressures would lead to higher 
aerodynamic heating and therefore to higher gradients. 


REFERENCES 


! Patterson, G. N., Molecular Flow of Gases, pp. 179-187; John Wiley & 
Sons, Inc., New York, 1956. 

? Schaaf, S. A., A Note on the Flat-Plate Drag Coefficient, University of 
California Inst. Eng. Res. Report HE-150-66, 1950. 

* Mirels, H., Estimate of Slip Effect on Compressible Laminar Buundary 
Layer Skin Friction, NACA Tech Note 2609, 1952. 


On the Effect of Air Pressure on 
Strouhal Number 


Anatol Roshko 


Guggenheim Aeronautical Laboratory, 
California Institute of Technology, Pasadena, Calif. 


August 25,1958 


5 em EXPERIMENTAL MEASUREMENTS of reference 1 show an 
effect of free-stream pressure on the Strouhal Number- 
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Reynolds Number relation for a vortex-shedding cylinder. Since 
it is unlikely that Buckingham’s z-Theorem needs re-examination, 
and since other possible parameters (Mach Number, etc.) would 
not appear to be relevant under the conditions of the experi- 
ments, one is at a loss to explain the results. It is interesting that 
if the Reynolds Numbers of reference 1 are divided by the pres- 
sures, normalizing on standard pressure, one obtains a single 
“standard” curve. Such a correction to Reynolds Number would 
be needed if account had not been taken of the density changes 
corresponding to the pressure changes. 
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Author’s Reply 


William D. Clifford 
National Aeronautics and Space Administration, 
Lewis Research Center, Cleveland, Ohio 


October 16, 1958 


IT THE EXPERIMENTAL PROGRAM presented on the effect of air 

pressure on Strouhal Number, the effect of pressure variation 
on density was taken into account in computing Reynolds Num- 
ber. Temperature was held constant at 80°F. Equipment was 
calibrated and extraneous effects including cylinder blockage and 
free-stream velocity profile were calculated. The data were not 
extensive but did indicate some interesting and temporarily unex- 
plainable trends. These data, along with other data showing 
similar temperature effects,! indicate that the Strouhal Number 
is not solely a function of Reynolds Number, as previously 


speculated. 
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A Generalized Porous-Wall 
‘“Couette-Type’’ Flow 


Kenneth R. Cramer 


Aeronautical Research Engineer, 
Aeronautical Research Laboratory, 


WADC, ARDC, Wright-Patterson Air Force Base, Ohio 
September 5, 1958 


oe it was observed that the two existing boundary- 
layer texts (references 1 and 2) did not contain a solution 
for the case of Couette flow with a constant, uniformly dis- 
tributed suction or blowing. Thus, the following analysis con- 
siders a ‘‘Couette-type’’ flow between a stationary flat surface and 
a slightly inclined flat plate moving at a constant velocity. In 
addition, the flow is subjected to a constant, uniformly distributed 
suction or blowing at the fixed surface. For this case, the 
boundary-layer Eqs. (1) and (2) are easily reduced to a simple 
ordinary differential equation in the following manner: 


u(du/dx) + v(du/dy) = U(dU/dx) + v(0*%u/dy?) (1) 
(0u/dx) + (dv/dy) = 0 (2) 
Boundary conditions: 


u=U 
y =0: u=0; v = const. = 


where b is negative for suction and positive for blowing. Letting 
a = u/U; 6 = the constant wall value = b/U 
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Velocity profiles. 


Fic. 1. 


= x/a; = ya; Re = Ua/sv and c = (a/v)(dU/dx) 
then Eqs. (1) and (2), neglecting the slight variation of a, become 
— + ¢ = (3) 
and = (4) 
with the boundary conditions: 
#@=0; 6=b/U =1 


Since Eq. (3) is independent of the x direction, the following 
ordinary differential equation results: 


— bRe(da/dj) +c = 0 (5) 
The solution of this simple equation is then Eq. (6): 
= Ky,e®Re3 + + Ke (6) 


Upon evaluating the integration constants, A; and Ao, with 
the boundary conditions, the final expression for the velocity 
profile, Eq. (7), is obtained: 
= {((c/#Re) — 1}/[1 — e®Re]} — 4 

(c/#Re)[¥ —1] +1 (7) 


i$ For the case of c = 0 and b = O, Eq. (7) reduces to the well 
Fe known linear velocity profile for plane Couette flow. Also, the 
linear velocity profile is obtained when c = 6Re, and this is the 
case where the streamwise velocity gradient is equal to the 
suction or blowing velocity. 

Fig. 1 presents several velocity profiles, and it can be seen that 
the combined effect of the velocity gradient and suction or blow- 
ing velocity can produce velocities greater than the plate velocity, 
or reverse flow. In general, adverse pressure gradients retard the 
a velocity profile in the vicinity of the fixed surface, and, for the 
particular pressure gradient shown, a value of Re = —1.594 
produces the separation velocity profile—i.e., the velocity gradi- 
ent at 9 = Ois zero. 

As Schlichting stated,! the analysis of Couette flow with a 
pressure gradient can be applied in the hydrodynamic theory of 
bearing lubrication. 
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Stability Criteria for Systems Having Complex- 
Coefficient Polynomial Characteristic 
Equations 


Louis G. Kaufman, || 


Project Scientist, Aeronautical Research Laboratory, 
ARDC, Wright-Patterson AFB, Ohio 
September 15, 1958 


WADC, 


gee woongpad IT IS NECESSARY to determine whether the roots 
of the characteristic equation of a system have negative or 
positive real parts, indicating, respectively, stability or instability 
of thesystem. Characteristic equations are generally of the form 


+ (A; + + (42 + +... 4 
A,+iB, =0 (1, 
/ 
where = V/-1 and A; and are real. 

Stability criteria have been developed because it is difficult 
to solve explicitly for the roots of Eq. (1) fer values of ” greater 
than two. The well-known Routh-Hurwitz criteria were de- 
veloped for the case of real coefficients (all B; = 0). Herbert 
von Bilharz developed stability criteria for the complex coefficient 
case in 1944;! Evelyn Frank independently arrived at the same 
results in 1946.2 The Bilharz-Frank criteria, for an nth degree 
polynomial, require the evaluation of m determinants, the largest 
of which is of order (2n — 1). 

The coefficients in Eq. (1) are frequently functions of some 
parameter which is allowed to vary. In these cases, it is neces- 
sary to perform many stability tests; this is quite laborious for 
n > 3 if the Bilharz-Frank criteria are used. In practice, it is 
usually found that » < 5. Hence, the following, considerably 
shorter, criteria were developed? for the cases of nm < 5. 

All of the roots of Eq. (1) will have negative real parts if and 
only if Di, De, ..., D, are all positive. The D,;’s are shown below 
for the case of m = 5. For n < 5 the expressions for the D,’s are 
identical to those shown below, but it must be noted that, for 
j > 1n, Dj; is not to be evaluated and A; = B; = 0. Thus, for 
n = 3, D; and D; are not to be evaluated and Ay = As = By = 
B; = Oare to be used in the evaluation of D,, Ds, and Ds. 

For A; #0: 


e 
dD, = Ai, = aBs cA, dD; = Ay Az —B, 
a b 
k mn 
da uv Ath w 
be ft and D; = 
Bo By A; @rs t 
where 
a = A,B, — B:, b = A,B; — B, 
c = A\A2 — A;, Ad; 
e = (aB,/A;) + Cy f (aB,/A,) om d 
g = (aA;3/A1) — h = A,B; — (aA;/A)) 
k= — (eA3/A)), m = —(eA;/A;) 
n = (eB./A;) g, p = (eBi/A\) —h 
qd A;B, A,B; 
A,B, “+ (cBs/A), r= A;B; A;B; (4;By/A1) 
s = d — (aB,/A)), t = A.A, + BoB, — A:A4 — 
= A,B.B;, + v = A,B,B; — dA; 


v= A,A3B; bA; 
If A; = 0, instability exists (by D2) unless Bz = 0 simultaneously- 
For A; = 0 = Bi: 


dD, => 0, Ds = 0, D; = rocls 


dD, A;(apA5By + horo goRo) doko Bs 


Ds = (noBsBs av? (doto foAsAs hor. AsBs 
Ropo)ao — + moro)bo + Az(moBs + qoAs)co + 
(goAsBs — qoro)di 
where 
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a = Bs — bo = ap AzB3; — — goAs — 
= rp AoBy + roBs + A;B; 


d = A3(acAs A;B;), + 
fi so( B;) ti Bs, go = A; AsA 
= 3.2 | B;By, ky = A 34 1; 


m, = AgAs + BsB;, mo = AzA5 — ho 
= + AsAsBs, ge = — A;Bs 
= A;?, = B,? 
= A;? 

The above stability criteria have proved to be quite useful to 
the author.* The criteria have been extensively checked alge- 
braically and numerically. 


REFERENCES 

1 Bilharz, Herbert, Bemerkung su einem Satze von Hurwilz, ZAMM, Band 
24, Heft 2, 1944 

2? Frank, Evelyn, On the Zeros of Polynomials with Complex Coe ficients, 
Bul. Amer. Math. Soc., Vol. 52, February, 1946 

3 Kaufman, Louis G., II, On the Stability of One-Dimensional Viscous, 
Compressible Flows: Sound Waves and Shock Waves, Ph.D. Thesis, Poly 
technic Institute of Brooklyn, June, 1958 


Experimental Investigations of the Effect of 
Heat Transfer on Hypersonic Turbulent 
Boundary-Layer Skin Friction* 


Eva M. Winkler and Moon H. Cha 
U.S. Naval Ordnance Laboratory, White Oak, Silver Spring, Md. 
September 15, 1958 


7 PRESENT TESTS were performed to investigate turbulent 
boundary-layer skin-friction and heat-transfer coefficients 
for various rates of heat transfer from the fluid to the wall and to 
provide a basis for analyzing the adequacy of existing theoretical 
approaches to the problem. The experimental studies were car- 
ried out on a cooled flat plate at a free-stream Mach Number of 
approximately 5. Tne free-stream conditions were selected such 


* The research discussed in this paper was sponsored jointly by the 
U.S. Navy Bureau of Ordnance and the U.S. Air Force. 
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Fic. 1. Variation of local skin-friction coefficient with momen- 
tum thickness Reynolds Number. 


that transition to naturally turbulent boundary layers occurred 
about 17 cm. from the leading edge of the plate. Downstream of 
the transition region, wall heat-transfer measurements and 
boundary-layer pressure and temperature surveys were obtained 
at several leading-edge distances and for three rates of heat trans- 
fer. The techniques employed and the data reduction procedure 
followed, in general, the methods described in reference 1; de- 
tails will be contained in a forthcoming NOL publication. 


TABLE 


M T,/T, Re, Reg 6* 
(Atm) (°K) x 107° (mm) (mm) 

5.21 7.11 383 0.861 2.72 2090 1.66 0.166 
5.14 7.08 345 0.875 3.36 2936 1.57 0.189 
5.20 7.21 366 0.840 4.07 3173 2.46 0.228 
5.26 7.31 364 0.845 5.04 3880 2.86 0.284 
5.29 7.09 363 0.845 5.94 4300 4.00 0.322 
4.98 5.42 383 0.757 2.29 1900 1.96 0.179 
5.18 7.01 399 0.744 2.58 1782 1.55 0.149 
5.20 7.05 374 0.753 3.81 2960 2.51 0.227 
5.24 6.99 379 0.773 4.55 3455 3.22 0.280 
5.24 6.55 384 0.765 5.11 3790 3.88 0.330 
5.17 6.77 455 0.613 2.01 1055 1.22 0.114 
5.16 6.71 476 0.586 2.57 1682 2.05 0.188 
5.10 6.84 468 0.578 2.73 1735 1.93 0.186 
5.11 7.20 496 0.565 3.27 2488 2.84 0.281 
5.20 7.51 467 0.585 3.57 2482 2.78 0.256 
5.12 7.00 465 0.605 4.22 3256 3.74 0.343 


~ 


~ 


Summary of Experimental Data and Comparison with Skin-Friction Coefficients Computed from Theoretical and Empirical Relations 


n c c c c c c c 
t t f t t f 
x10 x10 x10 x10 x «x 104 
q) (2) (3) (4) (5) (6) (7) 
9.9 14.65 14.87 13.93 16.32 12.70 10.20 15.30 
9.3 14.32 13.48 12.69 14.62 11.53 9.33 11.30 
8.2 13.46 13.43 12.08 13.32 10.95 9.08 10.71 
7.0 13,08 12.76 11.62 12.52 10.53 8.81 10.29 
9.7 13.35 13.48 15.17 17.86 14.01 11.31 13.77 
1.2 16.06 - 14.51 17.24 13.39 10.79 13.15 
9.9 12.48 11.65 13.55 14.92 12.25 10.10 11.98 
8.0 11.51 11.70 12.64 14.02 11.58 9.62 11.32 
6.7 oa 10.63 12.14 12.98 11.18 9.37 10.92 
6.9 14.70 15.47 15.90 21.60 15.02 12.18 14.74 
9.5 13.23 13.35 15.35 19.04 14.58 12.03 14.22 
9.8 12.36 12.79 14.84 17.74 14.19 11.95 13.80 
9.1 12.03 - ‘ 
8.7 10.54 10.75 13.94 16.42 13.14 11.16 12.78 


= momentum thickness 
= local skin-friction coefficient 


M = Mach number 
Pp = supply air pressure Experimental Data Theoretical Data 
if = supply air ‘temperature 
T = wall surface temperature (1) evaluated from boundary layer (3) Van Driest, ref. (2) 
Re = length Reynolds number velocity profile slope (4) Persh, ref. (3) 
4 (5) Colburn ) 
thickness Reysolds susber (2) evaluated from wall heat- (6) Nikuradse ) ref. (4) 
> = displacement thicksess transfer measurements (7) Schultz-Grunow) 
c 


5 


power-profile exponent 


Nore: Fluid properties involved in Rez, Reg, and cy are evaluated at the edge of the boundary layer. 
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The experimental results were compared with the theoretical 
predictions of Van Driest? and Persh’ and with a number of em- 
pirical formulas derived from incompressible data. To compute 
the skin-friction coefficients from the Colburn, Nikuradse, and 
Schultz-Grunow relations,‘ the properties were evaluated at the 
Rubesin and Johnson reference temperature (7’ method). © The 
results are given in Table 1 and illustrated in Fig. 1. Note that 
the theories have been plotted for constant values of 7,,/T> of 
0.60, 0.76, and 0.85 which correspond to the average values for 
the three groups of data. The following general conclusions may 
be drawn from the results: 

Good agreement (for most cases, within 3 per cent) was found 
between the skin-friction coefficients evaluated from the slopes of 
the boundary-layer velocity profiles and the values computed 
from the wall heat-transfer data, using Colburn’s version of the 
Reynolds analogy. 

The decrease in the value of the skin-friction coefficient with 
increasing momentum thickness Reynolds Number is consistent 
with the data reported in references 1 and 7. 

The decrease of the value of the skin friction with increasing 
momentum thickness Reynolds Number becomes more pro- 
nounced with increasing rate of heat transfer. 

The experimental data show that for a fixed value of the mo- 
mentum thickness Reynolds Number the skin-friction coefficient 
decreases with increasing rate of heat transfer to the flat plate. 
This behavior is opposite to that predicted by the theories.’ 

In comparison to the theories of references 2 and 3, the data 
show that the equations derived by Van Driest and Persh do not 
account correctly for the effect of heat transfer. The relations 
developed by Persh agree well for conditions of low heat transfer 
(T,,/T> ~ 0.85), but the predicted values become increasingly 
too large (up to 59 per cent) with increasing heat transfer. Van 
Driest’s prediction is too low for conditions of low heat transfer, 
is in fairly good agreement for moderate heat transfer (7/7) ~ 
0.76), and becomes too large (up to 34 per cent) for the high 
heat-transfer conditions. The same is true, more or less, for the 
interpolation formulas of Colburn, Nikuradse, and Schultz- 
Grunow. In general, the latter give the best fit to the present 
experimental data for the higher heat-transfer rates (7/7) ~0.6). 
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The Viscous Flow Near a Stagnation Point 
When the External Flow Has Uniform 
Vorticity 


J.T. Stuart 
National Physical Laboratory, Teddington, Middlesex, England 


October 2, 1958 


T° VIEW of the recent controversy between Li’ ® 4 and Glauert® 
on the nature of the solution of the boundary-layer equations 
when the external flow is rotational, it seems worthwhile to draw 


attention to a certain exact solution of the Navier-Stokes equa- 
tions which lends support to Glauert’s point of view. 

Consider the two-dimensional flow in the semi-infinite region 
bounded by the plane x = 0; let z be the distance normal to 
the surface; let « and w be the corresponding velocity com- 
ponents, and let » and p denote pressure and density, respectively. 
A solution of the inviscid equations of motion and of continuity 
is given by 

“= ax—ws, w= —az, 


p = —(1/2)pa%x? + 2) + const.) 


where a is a constant and w is the (constant) vorticity, namely, 
(Qw/dx) — (du/dz). 

It is known in the case of w = O that there is a related solution 
of the Navier-Stokes equations for which the flow outside the 
viscous boundary layer is of the form u = ax, w = —a(z — 6), 
where 6, is the displacement thickness. This solution shows that 
the effect of viscosity on the external flow is to displace the 
streamlines outward from the boundary by an amount 4. 
We wish to examine now if a solution of this kind exists when 
The Navier-Stokes and continuity 


(1) 


w is a nonzero constant. 
equations are 


u(du/dx) + w(du/dz) = —p~(dp/dx) + 
u(Ow/dx) + w(dw/ds) = —p-"(Op/ds) + (2) 
(Ou/ox) + (Ow/dz) = 0 \ 
where vy is the kinematic viscosity. If we now substitute 
u = axf'(n) — w(v/a)'”? g(n) 
w= —(av)'?f(n), = 2(a/v)'? (3) 
p = —(1/2)pa?x? — pav F(n) \ 
we find that f, g, and F must satisfy 
1 f" = 0 (4) 
g’ + fe’ —f'g =0 (5) 
F = (1/2)f? + f’ + const. (6) 


The boundary conditions are that the velocity shall be zero at 
n = O, and that the velocity shall tend to the solution (1) as 
n— © (aside from a possible displacement effect in the » co- 


ordinate). Thus, we have 
f(0) = f(0) = 0, = 1 (7) 
g(0) = 0, g(o) = 1 (8) 


It will be noted that formulas (1) and (3) may formally be 
generalized by replacing x by x + a; this, however, merely 
implies a shift of the stagnation point. 

Eq. (4) with the boundary conditions (7) defines the classical 
stagnation point flow, for which the function f(7) is known and 
tabulated.5 We now seek a solution of Eq. (5) subject to the 
conditions (8). As Glauert® and Rott? have pointed out in 
another connection, a particular solution of Eq. (5) is g = f” (7). 
The general solution of Eq. (5) therefore, is, 


g = exp fan + a] (9) 


where C, and C2 are constants. Glauert (Appendix )* has shown 


that 
\ 
(f ) exp; -f, fdn (an 1.335(n — 0.6479) (19) 


asn — ©. It follows, therefore, that the solution of Eq. (5) 
which satisfies the conditions (8) is given by C; = 0.7492, C. = 0, 
and is 


When » > ©, it can be seen from Eqs. (10) and (11) that 
— 0.6479, so that when s © 


u~ ax — w(z — 6), 
w~ —a(s — 6), (12) 
where 6; = 0.6479(v/a)'? is the displacement thickness. This 


result is similar to Eq. (1) but shows a displacement in the 
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coordinate. It, therefore, supports Glauert’s viewpoint’ that 
the proper condition at the edge of the boundary layer in the 
Blasius case is u — uo + w(y — 5), where y is the coordinate 
normal to the plate, and the oncoming stream is u = uo + wy. 

The skin friction on the wall can be calculated from Eq. (11) 
and is te = p(Ou/dz)o = pl1.2326 ax(a/v)"? — 0.6078]. 
The part of the shear which is proportional to w is seen to be 
numerically less at the wall than in the external stream.* 

It is of interest to note that the solution described above is a 
member of a class of solutions described by Lin.® 
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On Thermoelastic Stress-Strain Relations 
for Thin Isotropic Shells* 


P. M. Naghdi 

Professor of Engineering Science, University of California, 
Berkeley, Calif. 

September 12, 1958 


I THE ABSENCE of thermal effects, the formulation of suitable 
stress-strain relations in the linear theory of thin elastic shells 
(where the effects of transverse shear deformation and transverse 
normal stress are accounted for) has been recently carried out 
by E. Reissner! for axisymmetric deformation of shells of revolu- 
tion, and by the present author? for the general shell where the 
deformation is referred to the lines of curvature. It is also rele- 
vant to mention here that the usual formulation of problems of 
thermoelastic shells of isotropic materials—e.g., as in Melan and 
Parkus*—in the spirit of Love’s first approximation is defective 
in the sense that it does not conform to the requirement that the 
initially stress-free isotropic shell (in the absence of suitable edge 
constraints), should remain stress-free and undergo only a uniform 
dilatation, when subjected to a uniform temperature field. 

It is the purpose of the present note to extend the previous re- 
sults? and to modify the derived stress-strain relations so as to 
include the thermal effects in a manner consistent with the as- 
sumptions made in reference 2._ While there is no difficulty in 
carrying out the extension mentioned, in view of the considerable 
recent interest in thermoelastic problems of shells, it is desirable 
to record for future reference the necessary and relatively minor 
modifications in the stress-strain relations for isotropic materials 
which (on account of the assumed functional form of the trans- 
verse normal stress and the normal component of the displace- 
ment )? will be free from the defect mentioned above. 

For brevity we employ the same notation as in reference 2 and 
also recall the assumed form of the normal displacement W and 
transverse normal stress o¢ given by Eqs. (2.7)2f and (2.9)2—i.e., 


* The results presented here were obtained in the course of research spon- 
sored by the Air Force Office of Scientific Research under Contract AF 
18(603)—47, when the author was at the University of Michigan. 

+ The designation Eq. (2.7): refers to Eq. (2.7) of reference 2. 
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W =w + tw’ + (1) 


o¢ = [1 + (¢/R:)] [1 + X 
({ (3/2) (S/h) + (T/4h) [¢/(h/2))} 
[¢/(h/2)}?} + (1/2)qtH* {1 + (3/2) [¢/(h/2)} (2) 
(1/2) [¢/(h/2)]*} + 
(3/2) [¢/(h/2)] + (1/2) 
where the functions w’, w”, S, and T are specified by Eq. (3.5)s. 
In order to include the thermal effects, it is necessary to modify 
Eq. (3.1b)s, so as to read 
A = + o2 + +... (3) 
where the dots in Eq. (3) refer to the right-hand side of Eq. 
(3.1b)s, @ is the coefficient of the linear thermal expansion (which 
may depend on temperature and coordinates), and 6 = @(£;, £, ©) 
denotes the change in temperature from the initially stress-free 
temperature state. 
With the notation 


{awl = a 2) 
ai + { 


dé (4a) 


+ 
2)" ((1/2) [¢/(4/2)]5 — de 


(4c) 


and with A is given by Eq. (3), some of the various coefficients in 
Eq. (3.4). are modified as indicated below, 


6S: 6T: — Box 


and the first two of Eq. (3.5). will contain additional terms—i.e., 
t+... (6) 

= +. ..S 
but the last two of Eq. (3.5)2, as to be expected, remain unaltered. 
It follows that, in the presence of a temperature field, the modi- 


fied stress-strain relations (for normal components of stress re- 
sultants and stress couples as well as shear stress resultants) read 


= [An — +... 
— +... 
kr = [Ou — +... 
= [Oe — (6:¢/R2)] +. (7) 
(2/28) 
((1/a2) (0/dE)) 
= + 30[¢/(h/2)] Ox) +... 
An examination of Eqs. (4) and (7) reveals that when @ = const. 
(say, 6), then 


& 
| 


= —(3h/2) 


On = On = Or = = Ol 


8a) 
= (aO/Ri), Oo = (a6s/R2)\ 


and = esq aby, = ke = 0 (8b) 


It is also worth noting that, in order for the theory of Love's first 
approximation (in the presence of thermal effects) to be free of 
the defect mentioned earlier, it will suffice to assume for the 
normal displacement the form 


W =w (9) 


in which case, since 62x = 0, Eq. (7) simplifies considerably. 
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On Slip-Flow Heat Transfer to a Flat Plate 


Richard A. Oman and Richard A. Scheuing 


Research Department, f 
Grumman Aircraft Engineering Corp., Bethpage, NY. 


September 29,1958 


SYMBOLS 
Subscripts 

co = free stream 

oe = adiabatic stagnation 

aw = adiabatic wall 

= continuum.flow value 

fm = free-molecule-flow value 

i = incident to the surface 

R= re-emitted from the surface 

s = on the fluid side of the assumed slip discontinuity 
w = at or within the wall 


M“ POINTS OUT an expression for slip-flow heat-transfer 
rate which includes the shear work transmitted to the wall 
by the slipping fluid; viz., 


= k(OT/Oy), + wu, (Ou/dy), 


where & is the thermal conductivity and uw the viscosity. We do 
not yet know if the continuum expression for shear work will re- 
main valid, or whether new effects will appear in slip flow; but 
for the simplified analysis of this note we will assume that Eq. (1) 
is valid. 

From Maxwell's derivation of the slip velocity we take the 
approximate relationships: 


U, = (1;/2) + o)/2)u; (2) 
ui; = Us + (OU/Oy)s (3) 
= [(2 a)/a| A(Ou/Oy)s (4) 


where o is the Maxwell reflection coefficient, and \ the mean free 
path. 

We assume that the boundary-layer velocity profile is linear 
and that the effective boundary-layer thickness is the same as in 
continuum flow, being unaffected to first order by the change in 
the wall boundary condition. 

The slip-plane velocity gradient then becomes (0u/d0v), = 
(to — us)/6, where 6 is the effective boundary-layer thickness. 
We can find the slip velocity in terms of the Knudsen Number; 
Kn = 

= }[(2 — «)/o] Kn} /{1 + [(2 — «)/o] Kn} (5) 


Note that the validity of the above expressions becomes in- 
creasingly doubtful as An increases, and Kn>1 has no meaning. 


THE TEMPERATURE JUMP IN A MOvING FLurp* 
Kennard? derives a simple relationship between the temper- 
ature discontinuity at a boundary in a rarefied gas and the local 
temperature gradient in the gas. We modify this expression to 
account for the effect of motion in the gas. 

After Kennard, the energy incident on a wall from a stationary 
gas can be written: 


E; = + (1/2)R]T. + (1/2)R(0T/dy), 


where S is the mass flow rate of incident molecules. 

For a moving gas, the incident and re-emitted energies are in- 
creased by the kinetic energies of the incident and re-emitted 
molecules due to their mass motion. Thus: 

E; = S{ [ev + (1/2)R]Ts + (ui2/2)} + (1/2)k(0T/dy), 

Ew = S{lev + (1/2)R] Tw + (1 — o)? (ui?/2)} 

The net energy transfer from the fluid to the surface is the 
difference between the energies of the incident and re-emitted 
molecules; g = E; — Er. Using the accommodation coefficient 
a = (E; — Epr)/(E; — Ex), which we assume to be identical for 
* The authors are indebted to Dr. C. du P. Donaldson for pointing out 


that the conventional temperature-jump expression must be modified when 
used for a moving fluid. 


126 JOURNAL OF THE AERO/SPACE 


SCIENCES FEBRUARY 1959 


a 
(2) 
2 
= 
a © FREE MOLECULE THEORY 
T (REF. 5) 
a-o0=09 
Pr=0.72 
12 
fm 
2 
x 
10 
oO 
2 
& 9 +4 2 
x 
oO 
WwW > 
2 S 
] 
q W 
> 
(@) 2 4 6 8 1.0 
KNUDSEN NUMBER ~ \/8 
Fic. 1. Effect of Knudsen Number on heat-transfer rate and 


recovery factor. 


all degrees of freedom, we write 


q = a (S{ [ee + (1/2)R] (T. — Tx) + 
[(20 — o2)/2]u:2} + (1/2)R(OT/dy).) (6) 


Since g must also be equai to the energy crossing a plane which 
is outside the discontinuities at the surface, we can combine 
Eq. (1) with Eq. (6) in order to eliminate g. From the simple 
kinetic-theory relationship given by Kennard, 4» ~ 2X, and 
from Eqs. (2) and (4), we write the expression for the temperature 
jump in a moving rarefied gas: 


Ts — Tw = + 1)ep] {[2(2 — a)/a] X 
(07 /dv), + — a)/a] [o/(2 — o)| (7) 
EFFECT OF SLIP FLow ON HEAT-TRANSFER RATE 


Schlichting* gives the energy equation for continuum flow in a 
boundary layer in the form 
—(Op/0x) [c,(dT/du) + u] + cpp[(Pr — 1)/Pr] xX 
(dT /du) (0?u/dy?) = [(d?T/du®)k + p] (Qu/dy)? (8) 


A solution to this equation for Pr = 1 and (0p/dx) = 0 is 
d?T/du? = —yp/k. We assume an approximate solution of the 
form 


ad?T/du? = —(p/k) [1/f(Pr)] (9) 


where f(Pr) will be determined by comparison with the familiar 


continuum recovery-factor expression, 7 ~ V Pr. 
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We assume that Eq. (8)—and also the other boundary-layer 
equations—will remain valid well into the slip region, provided 
appropriate boundary conditions are applied to account for the 
effect of appreciable mean free path. The conclusions of several 
authors* ® support this assumption, especially for An < 0.2. 
We therefore apply the slip-flow boundary conditions T = T, at 
u= us,and 7 = at u = to the second integral of Eq. (9), 
and differentiate to get 
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(OT /Oy), = Ua(Ou/Oy)s (( Pr/2fcp) [1 — + 
{((T2 — Tw) + (Tw — — 
(10) 
We express 7), in terms of the independent variable » 
+ — To) (11) 


Incorporating Eq. (11) into Eq. (10), along with the temper- 
ature-jump relationship, (Eq. 7), we get 


Ua (Ou/Oy), [Pr/f(Pr)] (1 — (us/uo)|? —n — [8y/(¥ + 1) 1 — a), (2 — )* 
2c, 1 — }1 — [2y/(y + 1)] [((2 — @)/a] [o/(2 — (1/Pr)} 
Substituting Eq. (12) into Eq. (1), we find an expression for g, which becomes, when f is evaluated at g = 0, u./ua = 0, and ma, = 


Sir = te 


2VP, 


Oy u 


The ratio (q/g-) can easily be calculated from Eq. (13) because 
gc is g for = 0. Fig. 1 shows vs. Kn for several values 
of n, at Fr = 0.72, and with o = a = 0.9. The results show 
the expected behavior for small values of Kn. 
RECOVERY FACTOR 

An expression for 7 is obtained by setting g = 0 in Eq. (1) and 
using it in Eq. (12) to solve for may = 7: 
r= + 1)] [¢/(2 — + 

— + 2Pr(us/to) [1 — (14) 


An alternate approach toy is based on the approximation 


T,—To * Pr — us?) which is fairly accurate for no 
heat transfer, and for T,close to Tay, —i.e., small By defini- 
tion 

— Tao) = 


((Taw — Ts) + — 


r= [4y/(y + 1)] [o/(2 — o)] + 


V Pri — (15) 


The Replacement of Lifting Surfaces 
by Lifting Lines With Variable Position 


A. Van Tuyl 
U.S. Naval Ordnance Laboratory, White Oak, Silver Spring, Md. 
October 13, 1958 


I THE FIRST PART Of reference 1, C. Ferrari has investigated the 
replacement of a lifting surface in a uniform supersonic flow 
by a curved or straight lifting line. Linearized theory is used, 
and the lifting surface and lifting line are assumed to lie in the 
same plane and to have the same spanwise distribution of circu- 
lation. Let x be measured in the direction of the uniform flow, 
and let x, y, and z be such that the forward Mach cone from 
(x, y, 2) contains the lifting surface and the lifting line. The 
results of reference 1 show that when x is sufficiently large for 
the given y and z, there is a variable position of the lifting line 
which depends only on x such that the velocity potentials due to 
the lifting surface and lifting line, respectively, agree at (x, y, s 
up to and including terms of the order of 1/x*. A single lifting 
line can therefore be used for calculation of the downwash and 
sidewash to this order of accuracy at points of a given plane x = 
constant. 


1 — {1 — [2y¥/(y + 1)] — @)] — @)/al (1/Pr)} 


Eqs. (14) and (15) are identical for Pr = 1. Both forms show 
the expected behavior as increases toward the free- 
molecule value (0.55 for « = 0.9), increasing to a value greater 
than unity. The assumptions involved in the use of Eq. (1) 
therefore appear to account for the predominant slip-flow effects. 
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The purpose of the present note is to show that under condi- 
tions which are satisfied in certain cases of interest, the preceding 
results hold when 1/x4 is replaced by 1/x®. In particular, these 
conditions are satisfied when the lifting surface, lifting line, and 
loading are symmetrical with respect to the xs-plane. 

Because of the condition that the forward Mach cone from 
(x, y, z) contains the lifting surface and the lifting line, the pre- 
ceding results are seen to hold also when the flow at infinity is 
subsonic. The linearized flow calculations involved are the same 
in both the subsonic and supersonic cases. 

Let the flow at infinity be parallel to the x-axis, with velocity 
magnitude U,, and Mach Number 1/,. Let z be measured 
opposite to the direction of the downwash, and let a variable 
point in the xy-plane be denoted by (£, 7). Let S be a lifting 
surface with leading edge = and trailing edge = £(n), 
and let / be a curved or straight lifting line with equation = 
A+ &(n), —b <n» <b. The velocity potentials due to S and / 
and their trailing vortices will be denoted by ¢s and @), respec- 
tively. The circulation, I'(y), at a given span position is the same 
for Sand 1. 

When M,, > 1, and x, y, and z are such that the forward Mach 
cone from (x, y, z) contains S and /, we have 


| 
| 
+ 
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f° dn 
bs(x, 2) — (x, y, 2; r) f 
(y— LW (x 8? — — 0)? — Be? 


- — 7) dé (1) 
V(x — & — »)? — By — 9)? — 


[reference 1, Eq. (8) ], where B? = Ma* — 1, and y(é, 7) is the jump in ¢s across the xy-plane. Eq. (1) holds also when B? < 0. 


A = A(x, y, 5) be the solution of the equation 


s(x, ¥, 2) — Oi(x, y, 2; 4) = 0 


When x is sufficiently large with respect to y and z, we have 


A(x, 2) = + [m/(x — + [n/(x — +... 


b & b 
where = dn (& — &) f° I'(n)dn 


&2 
m = (3/2) — Aj)? — Eo? ] f° ronan 


and 


b \ b 
b / 


The condition for m to be independent of y and g is 


Je, a-h- 


New) (Ou/dE) dé = 0 (7) 


Given £(»), &(), and u(én), Eq. (7) is satisfied by an infinite number of functions (7). In particular, it is satisfied when S and / are 


symmetrical with respect to the xz-plane, and u(é, 7) is an even function of 7. 


Eqs. (3)—(7) reduce to equations for a straight lifting line 


when £() is set equal to zero. The quantity A is then equal to the abscissa of the center of pressure of S. When Eq. (7) is not satis- 
fied, it is seen that m depends only on y. We can then use / for the calculation of downwash, but not of sidewash, up to and including 


terms of the order of 1/x°. 


The downwash behind the straight lifting line & = 4, —b < n < b, ina supersonic flow is 
b b 
b 
| B4/[16m(x — d)*]} (y — 2) [(y — + (dP/dn)dn +... (8) 


when x, vy, and z are such that the forward Mach cone contains 
the lifting line. When z = 0, the first integral is interpreted as 
Cauchy’s principal value. This can be obtained, for example, 
from Eq. (21) of reference 2. When A in Eq. (8) is replaced by 
Eq. (3), we obtain an approximation to the downwash behind S 
which is accurate up to and including terms of the order of 1/x°. 

As an example, let S be a triangular flat plate at incidence in a 
supersonic flow, with = aln|, = 1,and B/a < 1. As in refer- 
ence 3, Eq. (24), we have 


u(é, = CV — a*,? (9) 


where C depends only on M,,, U,,, and the angle of incidence. 
The values of downwash obtained from Eqs. (3), (8), and (9) 
agree closely with the exact values given by linearized theory at 
points for which the forward Mach cone contains the entire lift- 
ing surface. The values of \ obtained from Eqs. (3) and (9) lie 
between 2/3 and an upper value which depends on B/a for x 
larger than 1 + B/a. The upper value remains between 0.777 and 
0.786 as B/a varies between 0 and 1. This is consistent with a 
result of Mirels and Haefeli,? who found that good agreement 
between a triangular flat plate at incidence and a straight lifting 
line is obtained for moderate values of x when the lifting line 
is at 3/, the distance from the vertex to the trailing edge. 

It is pointed out in reference 1 that the best straight lifting line 
lies between the center of pressure of S and the trailing edge when 
x is sufficiently large, and that this is in agreement with the pre- 
ceding result of Mirels and Haefeli. An explanation of Mirels 
and Haefeli’s result given by D. C. M. Leslie (reference 4, p. 
299) is seen not to be correct, since the fixed straight lifting line 
which gives the downwash behind S up to and including terms 
of the order of 1/x* has been shown (reference 1) to pass 
through the center of pressure. 
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be double or triple spaced on one side of white paper sheets, 
consecutively numbered. There should be wide margins to allow 
for the marking of directions to the printer. Correcting, chang- 
ing, or adding to papers after they are in type is costly. It is, 
therefore, imperative that papers submitted be in final form. 
Typographical errors may be corrected on proofs, but if authors 
wish to add material, they may do so at their own expense. In 
mailing, drawings may be rolled, but manuscripts should be sent 
flat. Send by first-class mail (register if you wish for your own 
protection) to the Editorial Office, Institute of the Aeronautical 
Sciences, 2 East 64th Street, New York 21, N.Y. All manuscripts 
will be examined by the Editorial Committee and by the Editor. 
# Authors will be advised as promptly as possible whether the paper 
is acceptable for publication. 

TitLES: The title of the paper should be brief. The name 
and initials of the author should be written as he prefers. The 
use of the full name of an author is advocated because of the fre- 
quent duplication of initials and surnames which sometimes makes 
it difficult to establish the identity of the author. This is par- 
ticularly important for large annual indexing and abstracting 
services. All titles and degrees or honors are omitted. The name 
of the organization with which the author is associated should 
be placed after his name on a separate line. The date on which 
the paper is received will be inserted by the Editor. The author's 
title or position should be indicated in a footnote. 

SUMMARIES OR AsBsTRACTS: An abstract to be printed at the 
beginning should accompany each article. It should be adequate 
as an index and asa summary. It should contain a statement of 
major conclusions reached, since summaries in many cases con- 
stitute the only source of information used in compiling scien- 
tific reference indexes. Abstracts printed in other journals, espe- 
cially foreign, in most cases consist of summaries from printed 
papers. The summary should explain as adequately as possible 
the major conclusions to a nonspecialist in the subject and should 
contain from 100 to 300 words, depending on the length of the 
paper. 

SUBHEADINGS: Subheadings should be inserted by the author 
at frequent intervals. The work of editorial preparation will be 
simplified by the author’s provision of many subheadings. 

MatrerR UsuaL_y DeLetep: Photographs or illustrations of 
little technical interest and not showing advances in general 
practice. Too detailed tabular matter (general results of such 
tables may be included in the text). Lengthy descriptions of 
materials or processes or of preliminary experiments or theories 
that preceded final results; salient features only are of interest. 

REFERENCES: References should be numbered consecutively 
and grouped together at the end of the manuscript, with only 
the corresponding number being mentioned in the text. The 
arrangernent should be as follows. For books: ' Durand, W. F., 
Aerodynamic Theory, 1st Ed., Vol. 1, p. 23; Julius Springer, 
Berlin, 1934. For magazines: 'England, Crawford, 
A. B., and Mumford, W. W., Some Results of a Study of Ultra- 
Short-Wave Transmission Phenomenon, Proc. IRE, Vol. 20, 
No. 12, pp. 481-482, March, 1933. Please give author, title, 
edition, volume, number, page, publisher, and place and date of 
publication as indicated. Omission of one required fact causes 
much extra editorial work and possible inaccuracies. 


ILLUSTRATIONS: Illustrations should accompany manuscripts, 
and each should always be referred to in the text by number. 
Drawings or graphs should not be larger than 12 X 16 inches, 
and must be made with jet black India ink on white paper or 
tracing cloth, the latter being preferred. Do not use typewriter 
for lettering. The smallest lettering on 8 X 10 inch figures should 
be no less than '/, inch high. Cross-section paper (white with 
black lines) may be used, but it should not have more than 4 lines 
per inch. If finer ruled paper is used, the major division lines 
should be drawn in with black ink, omitting the finer divisicns. 
In the case of finely ruled paper, only blue-lined paper can be 
accepted. Tracing paper and blueprints are not acceptable. 
Lettering and all markings must be large enough to be readable 
after reduction to single-column width (35/,.in.). Mail rolled or 
flat; never fold. Drawings that cannot be reproduced (including 
pencil drawings) will be returned to the author for redrawing, 
thus delaying publication of the paper. Photographs should be 
distinct and show clear black and white contrasts. They must 
be on glossy white paper. Avoid round and oval photographs. 


CAPTIONS AND LEGENDS: Legends or captions must accompany 
each drawing or photograph submitted. If written on the draw- 
ing or photograph, they should be placed below and well outside 
the part to be reproduced. Each table should have a caption 
such as Table 1, Table 2, Table 3, etc. Captions should be com- 
plete in themselves so as to make the data intelligible to the 
reader without reference to the text. A duplicate list of captions 
for figures should be included as the last page of the manuscript. 
Use ‘‘Fig. 1’’ (not Figure 1), ‘Figs. 3 and 4,” etc., in both the text 
and the numbering of illustrations. In the text, ‘Eq. (1)” or 
“Eqs. (1) and (2)’’ is used, not ‘‘Equation (1)" In captions, 
legends, and in table headings, write all words in full; do not 
abbreviate, except for ‘‘Fig.”’ and ‘‘Eq.” 


MATHEMATICAL WorK: Formulas may be typewritten or 
carefully written in pen and ink, the writing to be large enough so 
that it can be marked for the printer. Considerable space for mark- 
ing should be allowed above and below all equations. All compli- 
cated equations should be repeated on separate sheets with plenty 
of space left for marking. The solidus should be used for sim- 
ple fractions appearing within the text. Make all expressions 
clear to the typesetter. Greek letters used in formulas should 
be clearly designated by name in the margin of the manu- 
script. The difference between capital and lower-case letters 
should be clearly distinguished and care taken to avoid confu- 
sion between zero (0) and the letter (0), between the numeral 
(one) and the letter (ell) and the prime (’), between alpha and 
a, kappa and k, u and mu, v and nu, n and eta. All sub- 
scripts and exponents should be clearly distinguished. Avoid 
complicated exponents and subscripts. Dots and bars over 
letters or mathematical expressions should also be avoided. 
When it is necessary to repeat a complicated expression, it should 
be represented by some convenient symbol. 


SyMBOLS AND ABBREVIATIONS: The symbols recommended in 
the American Standard Association ‘Letter Symbols for Aero- 
nautical Sciences,’"” ASA Y10.7—1954, should be used wherever 
practicable. All symbols should be clearly written and carefully 
checked. Standard abbreviations should be used, and it should 
be noted that most abbreviations are lower case, such as m.p.h., 
b.m.e.p., i.b.p., b.hp., hp., ete. 
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